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ABSTRACT 



Three methods for weighting an exponential regression 
model to estimate discrete reliability growth were derived 
and tested. The first method systematically applies greater 
weight to test phases whose estimates have less variability. 
The second method similarly applies heavier weight to the 
most recent test phase estimate. The third method allows 
the user to choose the weighting scheme. 

These methods were evaluated against eight patterns of 
actual reliability by altering a previously developed Monte- 
Carlo simulation. Their performance was then compared to 
the unweighted exponential regression and Maximum Likelihood 
Estimate With Discounting (MLEWD) models. The second 
weighting method appears to perform the best under rather 
general constraints. Also, comparison is made of the least 
squares estimates for reliability growth using two different 
unbiased estimates for the negative of the natural logarithm 
of the failure rate. 
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INTRODUCTION 



Systems and equipment often go through several phases of 
development during which the reliability is assessed and the 
components or design are appropriately modified. From these 
early assessments, project managers need to make predictions 
regarding the ability of the system or equipment to meet 
reliability specifications by the contracted acceptance or 
delivery date. Some conventional estimators of reliability 
disregard data related to immature systems and use only that 
data obtained from the final product in verifying 
reliability. This practice is extremely inefficient, costly 
and occasionally, infeasible. Using all available test data 
would not only permit a manager to monitor the development 
process more successfully but at a substantial savings to 
all parties involved. 

Reliability growth models are useful in estimating the 
reliability of both immature and mature systems. Continuous 
and discrete versions of models exist. Models which employ 
attribute test data are discrete while those based on time 
to failure data are continuous. Only discrete models will 
be addressed in this thesis. Both of these versions make 
use of all available test data and are often used in 
conjunction with a technique known as failure discounting. 
This is a process whereby an increasing fraction of a 
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failure is repeatedly diminished as greater amounts of 
subsequent data is accumulated without repeat of the same 
failure . 

Models evaluating data which have been adjusted using a 
failure discounting scheme have performed well in estimating 
actual reliability [Ref. l:pp. 52-55]. Relying on 
discounting techniques has its drawbacks however. In order 
to employ the straight percent discounting method, for 
example, the user must select two parameters, N and F. F is 
the fraction each failure will be discounted and N is the 
number of successive successes which must be achieved prior 
to applying the discounting method. No rules have been 
developed as to how these parameters should be selected and 
no guidance other than to "use good engineering judgment" 
[Ref. l:p. 50] has been offered on the subject. 

The objective of this paper is to provide the user of a 
reliability growth model with an alternative to failure 
discounting which will produce comparable results. To this 
end, the Monte-Carlo simulation developed by Captain James 
Drake and modified by Captain James Chandler was altered to 
include a weighted exponential regression model . This model 
can use one of three methods of weighting. 

Method one is based on the idea that the variances of 
the observables are not all equal and hence some are more 
reliable than others. To accommodate this fact, weights in 
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this method are "heavier" for observables with lower 
variances. Method two systematically gives more weight to 
the most recent test phase motivated by the idea that the 
most recent version of the system is likely to be the most 
reliable. Method three is an extension of this idea. In 
this method, the user is permitted to subjectively select 
how much weight he desires to give each phase. This might 
produce successful results if the user has an intimate 
knowledge of the development process. He could then 
possibly choose to give little weight, say, to a phase he 
knew was conducted poorly or contained problems that he was 
certain were since corrected. Each of these methods was 
evaluated and compared to previously developed models using 
eight different actual reliability growth patterns. These 
comparisons are contained in Chapter VI . 

As a further modification to previous work which 
primarily explored the "test-fix-test" methodology, i.e. 
only one failure per phase of testing permitted, this thesis 
investigates the effect of using an alternative unbiased 
estimator of the exponential regression parameter of a 
"test-find-test" scenario. 

The following chapter will address previous work in the 
area of reliability growth models. The motivation for 
altering the exponential regression parameter is discussed 
in Chapter III and the weighted Regression Model and the 
derivation of the three methods of weighting are developed 
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in Chapter IV. A brief description of the models with 
which the weighted model will be compared is presented in 
Chapter V. The comparisons under varying patterns of actual 
growth is treated in Chapter VI. A summary, conclusions and 
recommendations for further study are discussed in Chapter 
VII. 
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II. PREVIOUS WORK 



A . BACKGROUND 

This thesis is the third in a recent sequence on the 
subject of reliability growth. The initial work was 
completed by Captain James Drake [Ref. 1] . In his study. 
Captain Drake developed a Fortran program to evaluate the 
performance of three reliability growth models; namely, the 
maximum likelihood estimate model, the exponential 
regression model, and the weighted average model. These 
models were evaluated in conjunction with two failure 
discounting methods, the confidence limit (C.L.) method and 
the straight percent discounting method. The C.L. method is 
referred to as the Lloyd discounting method in Drake's 
thesis. In addition. Captain Drake varied the parameters of 
each of the discounting methods to evaluate their effect on 
the performance of these models. The "actual" or known 
reliability growth pattern against which these models were 
compared was generated by Monte-Carlo simulation. In this 
simulation, the user inputs the reliability during the first 
phase of testing of each component which may possibly fail 
(i.e., of each potential failure cause). The program then 
generates an "actual" growth pattern from these values for 
the remaining phases of testing. The user of this program 
relinquishes strict control of the actual reliability growth 
pattern beyond the first phase and, although manipulation of 
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another parameter (specifically, "FRIMP" ) can produce the 
pattern in a general desired form, the simulation is 
incapable of producing a decreasing pattern. 

Captain James Chandler modified Captain Drake's program 
to permit the user complete control of the actual 
reliability growth pattern at each phase of development. In 
his version of the program, the user is required to input 
the reliability of each potential failure cause at each 
phase. Additionally, if more than one failure is permitted 
per phase of testing, the value must be replicated for each 
failure in that phase [Ref. 2] . To reduce the amount of 
required user input, this process was modified in this 
thesis so that the amount of information required per phase 
is the same regardless of the number of failures per phase 
permitted. The required information is contained in the 
User's Guide, Appendix A to this thesis. 

Captain Chandler also altered the C.L. failure 
discounting method to allow the user control of the discount 
interval. This was motivated by the results of the original 
thesis which indicated that employment of the C.L. method 
resulted in overly optimistic estimates regardless of the 
model used or the actual reliability pattern generated. 
This was attributed to the fact that the original C.L. 
method diminishes the weight of a failure after each 
following success. 
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B. RESULTS OF PREVIOUS WORK 



The weighted average reliability growth model was 
eliminated as a potential model in the early stages of the 
original analysis. This was due to its consistent 
overestimation of actual reliability which was only 
amplified by employment of either of the failure discounting 
methods. This model was discarded in subsequent work as 
well . 

The remaining two models accurately tracked a wide range 
of reliability growth patterns. The Maximum Likelihood 
Estimate with Discounting model, (MLEWD), generally "tended 
to underestimate actual reliability in early phases and 
slowly converge to the actual value with increased test 
data" [Ref. l;p. 48]. This model also exhibited smallest 
variance of all models evaluated. The negative aspect of 
the MLEWD model is that the choice of discounting parameters 
is critical to successful implementation and guidance in 
their selection is non-existent. The potentially drastic 
effects of different parameter choices is discussed in 
Chapter V. 

The exponential regression model also performed well 
against various actual reliability growth patterns and was 
far less sensitive to the choice of discounting parameters 
than the MLEWD model. The difficulty with this model is 
that it is highly variable in the early phases of testing 
although it generally stabilized after four phases. This 
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was considered acceptable since most development processes 
allow for more extensive testing. 

The C.L. discounting method was not recommended for 
employment with either model as its use only occasionally 
produced results comparable to the straight percent 
discounting method without possessing its corresponding 
flexibility. 

Captain Chandler evaluated these two models along with 
the standard single phase maximum likelihood estimate 
against eight reliability patterns. His results indicated 
that both the MLEWD and the exponential regression model 
were superior in all respects to the standard estimate. 
Additionally he found that the modification made to the C.L. 
discounting method resulted in instances of superior 
performance to the original C.L. method (although not to the 
straight percent discounting method). The remainder of his 
work substantiated the conclusions drawn by Captain Drake 
and included additional observations regarding the case of 
declining reliability. 

Both authors conducted all simulations for the test-fix- 
test scenario. 
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III. MOTIVATION FOR IMPLEMENTING Y *. 



The exponential regression reliability growth model was 
developed by H. Chernoff and W. M. Woods. Its derivation 
is fully detailed in Ref. 3 and Chapter IV. In this 
model the reliability after the change is modeled by 

Rj^=l-e“ . For the purposes of this chapter it is 

sufficient to understand that the model estimates the 
reliability in Phase k as; 






(3.1) 



where Rj^ is the reliability estimate, and the estimates 

A ... 4-h 

and for a and 6 at the conclusion of testing in the k 

phase are obtained using linear regression methods and an 
unbiased estimator for ( • 

The unbiased estimator is; 



_ n + l/2+l/3 + ...l/(.Yy^-l) 

- jo 



for Xji,>2 
Otherwise 



(3.2) 



where Xjj^ is the number of trials between the (j-1)®^ 
failure and the failure (including the failure) in 

the k^^ phase. 



constant throughout 

•cnis thesis. B is used for 6 within equations. 
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Let Fj^ = the number of failures in the phase. Since 

^l,k' ^2,k'***^F ,k • independent random variables then 

}c 

Yk = {Yxk^Yik^-Yf;)\h ( 3 . 3 ) 

is also unbiased. 

Previous work incorporates this estimator of (a + 
however an unbiased estimator which has minimum variance of 
all estimators has since been developed by W. M. Woods; 
this estimator is as follows: 

y/ = \jF, + \1{F, +1) + - 1) 

Replacing Yj^ with y]^ in the exponential regression model 
provides a more accurate estimate of reliability. 

Table 1 is the result of performing ten thousand 
replications of one phase negative binomial (Pascal) test 
data with actual reliability R and comparing the model to 
this value of R employing Yj^ and again using Yj^. This table 
demonstrates that while reliability predictions using 
y]^ and Yj. are both conservative, those using more closely 
estimate actual reliability, have smaller variances and a 
smaller mean square error. 

The parameters used in this simulation are as follows; 

Fj, = The number of failures in phase k was fixed at 3 for 
all cases 

Xp = The number of trials to the Fj^ failure. 
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N = The number of replications = 10,000 

k = The number of phases = 1 

and are as in equations 3.3 and 3.4 respectively. 

TABLE 1. COMPARISON OF Yj^ VERSUS Yj^ AS THE ESTIMATOR OF Aj^ 



Actual 

Reliability 


A 

R, 


A. 

R, 


4 


2 




1 ^ 2 


.5 


.430058 


.452505 


.053019 


.051689 


.053714 


.053798 


.6 


.548(M5 


.551300 


.049314 


.047326 


.031683 


.049269 


.8 


.737407 


,762389 


.026653 


.024093 


.029267 


.024800 


.9 


.8'2237 


.876379 


.011429 


)( I9N99 


.012034 


.009436 


.95 


.933223 


.936339 


.0(i48S9 


.(KM 163 


.(KM 280 


.(K)3227 


.99 


,')S5l 1 1 


9S6n7l 


.001414 


.001348 


.002313 


.001747 



Though these results prove that Yj^ is superior to Yj^ as 
an estimator of the coefficient of the exponential term, the 
impact of this improvement is disappointingly trivial. 
Figures la - Ih depict the results of employing Y* and Y to 
the exponential regression model for actual reliability 
patterns 1-8. As one can see from both the Table and the 
figures, the effect of this change takes place in the third 
decimal value. In practice, changes of this magnitude are 
hardly useful. Still, since y]^ is an improvement over 
Yj^ it is retained as the method of estimating (a and 6k) 
throughout the remainder of this thesis. 

Note: In all figures, A=actual; *=Y* and - = Y 
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RELIABILITY 




Figure lA, Y* vs. Y, Pattern 1 
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RELIABILITY 




Figure 1C. Y* vs. Y, Pattern 3 




Figure ID. Y* vs. Y, Pattern 4 
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RELIABILITY RELIABILITY 

0.875 0.880 0.885 0.890 0.895 _0^4_^0.5 0.6 0.7 



I 




4 6 8 

PHASE (3 FAILURES PER PHASE) 

Figure IE. Y* vs. Y, Pattern 5 
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PHASE (3 FAILURES PER PHASE) 

Figure IF. Y* vs. Y, Pattern 6 
14 
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RELIABILITY RELIABILITY 

0.36 0.37 0.38 0.39 0.40 0.55 0.56 0.57 0.58 0.59 0.60 




Figure IG. Y* vs. Y, Pattern 7 




Figure IH. Y* vs. Y, Pattern 8 
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IV. THE WEIGHTED EXPONENTIAL REGRESSION MODEL 



A. THE EXPONENTIAL REGRESSION MODEL 

The exponential regression model obtains sequentially 
updated estimates Rj^ of the reliability Rj^ which denotes 
true reliability after the phase. The basic model for 

Rj, is: 

1 - (4.1) 

In the exponential regression model, linear regression is 
used to estimate Aj^s (ct+6k) by Aj. where 

The exponential regression estimate of reliability is then: 



Ru = I — e 



- (** + 



A =1,2,. 






(4.3) 



and 2 + 6j^k is estimated sequentially at each phase. This 

results in a model which is capable of tracking changing 
reliability. In order to estimate the parameters and 
one must first calculate the unbiased estimator, As 

discussed in Chapter III, is the unbiased estimator with 
minimum variance and has been implemented in this version of 
the model. 

Let; 

Fj^ = The number of failures in the k^^ phase 

Xp = The total number of trials to and including the 
k failure in the k^“ phase. 
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Then: 




As an example, if the testing in the fifth phase 
continues until three failures occur and the third failure 
occurs on the eighth trial then; 

Y 5 = 1/3 + 1/4 + 1/5 + 1/6 + 1/7 = 1.092857 

The least squares estimates for and a ^ are ; 



Replacing the unknown parameters, their 



k 




(4.5) 




(4.6) 



where k = (I + 2 + 3 + ...k)lk 

L = Y=\ + ...Y\)lk 



estimates and in equation (4.3) yields; 




for A > .1 
for k=\ 



(4.7) 
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B. THE WEIGHTED EXPONENTIAL REGRESSION MODEL 

The weighted exponential regression model is identical 
to the unweighted exponential regression model above with 
the following exceptions: 

Let: 



Wj^ = the weight applied to the estimate in phase i. 
Replace all occurrences of k in the calculation of the 

unweighted model with: 

k 

- t) 



k 



W 




(4.8) 



i=\ 

and replace Yj^ everywhere with: 

k 




/=! 

Making these substitutions, the estimates of 
now: 



(4.9) 

a and B are 



^[(/-/:Jx r* X wj 

^ (4-10) 

y (/ - kj- X Wi 
/=! 



and the estimate of reliability is: 



(4.11) 



R, = I • ~ ^ for k> I 

*'*' 1^1— e ‘"=1 — e * for k = \ as before 



(4.12) 
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C. EXAMPLE OF THE WEIGHTED EXPONENTIAL REGRESSION MODEL 

The following example is offered as a means of 
clarifying the application of this model. 

Suppose that a system undergoes ten phases of testing 
and a phase is terminated upon occurrence of the third 
failure in each phase (i.e., F2=3, F2=3.... F2 q= 3) . Suppose 
further that the weights given to each phase are as follows 
(See Table lA) : 

TABLE lA. EXAMPLE OF WEIGHT DATA 



Phase 


1 


2 


3 


4 


5 


6 


7 


8 


9 


10 


Weight 


.03 


.03 


.03 


.03 


.03 


.15 


.15 


.15 


.20 


.20 



These weights have been arbitrarily selected for the 
purposes of this example. Actual calculation and/or 
selection of weights will be discussed in section D of this 
chapter. 

Tables 2 and 2A represent the test data collected on 
this fictitious system and the resulting calculations based 
on this data: 
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TABLE 2 



EXAMPLE PARAMETER CALCULATIONS BY PHASE 



Pliase(k) 


fi Trials to 


y\ 


A. 








I 


4 




1 


.333333 


S-. 


* 


2 


4 


.333333 


1.5 


.333333 


.124 X in ^ 


.33333 1 


3 


7 


.949995 


7 


.538SSS 


.30S333 


-.077777 


4 


10 


l.32S9(>7 


2.5 


.73(>4oS 


. 3o^ >3.>o 


-.104484 


5 


4 


.333333 


3 


.055793 


.0995(>3 


.357104 


6 


9 


l.217S.'6 


10.^ 


2.24S4I1 


-.0933450 


3.22S533 


7 


4 


.333333 


■).333 


1.0lon5n 


•.|(>S3.t4 


2.021357 


8 


35 


2.6IS2n3 


9.0 


l.S02t»So 


-.125571 


2.99223(» 


9 


112 


3.791201 


1 1.25 


2.SO')S87 


-.123(.SS 


4. 20 1384 


10 


170 


4.2447S5 


I I.o 


3.(»9(>S()(’ 


-.13.^9. IS 


1 4..V123nS 



'' A '' 

and reliabliliiy estimates, /?*=!— 



TABLE 2A. EXAMPLE RELIABILITY ESTIMATES 



Pliase (k) 


1 


1 


3 


4 


5 


0 


7 


8 


9 


10 


R. 


.36840‘ 


.'37925' 


.43753^ 


.50047( 


.504331 


.92104V 


.89404' 


.891S4^ 


.95503: 


.96368J 


R. 


.3984 


.42SIO' 


.4S07‘r 


.53924' 


.(i0')94> 


.70208^ 


.79812-^ 


.89996' 


.95009( 


.990rM( 



* = not computed in the first phase. See equation (4.7). 
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D. METHODS OF WEIGHTING 



1 . Method One . 

In the unweighted exponential regression model the 
estimate of 6 is the Gauss-Markov estimate when the variance 
in each phase is equal. When the variance is not equal from 
phase to phase ”the variance of B is unnecessarily large.” 
[Ref. 5] . This variance in parameter 6 accounts for the 
demonstrated higher reliability estimate variance evidenced 
in previous work [Ref. 1]. In order to diminish this 
variability each phase is systematically accorded a weight 
determined by the estimate of the variance of the phase: 



= the weight allocated in phase k 




(4.13) 



using weights derived in this manner, in equation (4.9) 

becomes: 

k 

k (4. 14) 

i=i 



which is unbiased and has minimum variance among all 
unbiased estimators of B [Ref. 5]. 

Figure 2 demonstrates the effect of weighting the 
exponential regression model by this method. 
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As one can see, the weighted model more closely 
estimates actual reliability, particularly after the third 
phase of testing. Chapter VI contains more detailed 
comparisons of models. 



A=ACTUAL, R=REGRESSION, 1=METH.1 



o 




Figure 2. The Weighted Exponential Regression Model 
in Conjunction with Method One Weights. 
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2 . Method Two 



The second method of weighting applied to the 
exponential regression model was less theoretically derived. 
This method systematically allots more weight to the most 
recent phase. The rationale for the implementation of such 
a method is that since the aim of the development process is 
to improve the system, the reliability in the latter phases 
will be greater than in earlier phases. While this is not 
necessarily the case, the method none the less has intuitive 
appeal . 

The weights used in this second method are 
calculated for each phase as follows: 




1=1 



( 4 . 15 ) 



To see why this expression creates weights which 
increase with phase, recall that an estimate of the variance 
of the negative binomial distributed data is; 



a2 



= an estimate of variance in phase k = 



jFk) X iPk) 



where : 

Fj^ = # of failures in phase k, as before, and 

A 

P]^ = The maximum likelihood estimate of the 

reliability in phase k. 



( 4 . 16 ) 
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Assume, as in our example in section C, that the 
number of failures in each phase is constant, then the term 
F]^ may be ignored. One can readily see that as pj^ 
increases, increases and therefore, Wj. increases. 

Figure 3 demonstrates the effect of weighting the 
exponential regression model by this method. The model 
weighted in this manner more accurately tracks actual 
reliability than both the unweighted and method one weighted 
model. More extensive comparisons are drawn in chapter VI. 



A=ACTUAL. R=REGRESSION, 1 =METH.1 ,2=METH.2 




Figure 3. The Weighted Exponential Regression Model 
in Conjunction with Method Two Weights 
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3 . Method Three 



The third and final method of weighting the 
exponential regression model has no roots in theoretical 
mathematics. This "method" is simply to let the user select 
the amount of weight he desires to assign to each phase. 
His only constraint is that the sum of these weights over 
all phases equals unity. Selecting the weights for this 
method is currently analogous to specifying the parameters 
for use in failure discounting. Unless one possesses 
extensive knowledge and intimate familiarity with the system 
being evaluated as well as its developmental history, use of 
this method is not advisable. 

To demonstrate the radical behavior of the model 
when this method of weighting is employed, five different 
cases of weights were evaluated. The weights for each case 
are listed in Table 3. 

Figure 4 is an example of the behavior of the model 
when the weights listed as case 3 are employed. As one can 
see, this selection of weights resulted in an extremely poor 
performance of the model against the pattern of actual 
reliability. As with the previous two methods, further 
comparisons are contained in Chapter VI. 
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TABLE 3. SAMPLE METHOD THREE WEIGHTS 



Phase(k) 


case 1 


case 2 


case 3 


case 4 


case 5, 


1 


.014285714 


.014285714 


.03 


.333333333 


.0025 


2 


.014285714 


.014285714 


.03 


.333333333 


.0025 


3 


.014285714 


.014285714 


.03 


.333333333 


.005 


4 


.014285714 


.014285714 


.03 


.333333333 


.0025 


5 


.014285714 


.014285714 


.03 


.333333333 


.01 


6 


.014285714 


.014285714 


.15 


.333333333 


.025 


7 


.014285714 


.014285714 


.15 


.2 


.05 


8 


.2 


'y 

.J 


.15 


.2 


.1 


9 


.3 


.3 


.2 


.2 


y 

o 


10 


.4 


■y 

.J 


.2 


.2 


.5 



A=ACTUAL, R=REGRESSION. 1 =METH.1 ,2=METH.2,3=METH.3 




PHASE (1 FAILURE PER PHASE) 

Figure 4. The Weighted Exponential Regression Model 
in Conjunction with Method Three Weights. 
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As a comparison, the weights for each of the methods 
depicted in Figure 4 are listed in Table 4. 

TABLE 4. COMPARISON OF METHODS ONE, TWO AND THREE WEIGHTS 
OVER A TEN PHASE SIMULATION. 



riiase 


1 


2 


3 


4 


5 


6 


7 


8 


9 


lU 


Method 1 
Weight.s 


I.O 


•5(M) 


.166 


.07583 


.»i8()63 


.06247 


.06827 


.01891 


.00318 


.00982 


Method 2 
Weights 


I.O 


.500 


.555 


.5,3017', 


.207695 


.279013 


.192177 


.4()5(U3 


.705975 


.695464 


Method 3 
Weights 




.03 


.03 


,03 


.03 


.15 


.15 


,I.S 


2 


.2 



The values listed are the amount of weight given to 
that particular phase when it is the current phase. Thus, 
for example, in Phase 9, method 1 allots .0031794 to the 9th 
phase and 1- .0031794 = .9968206 to the previous eight 
phases. Method 2 allots a much greater weight to phase 9, 
i.e., .705975, and only 1-. 705975 = .294025 to the previous 
eight phases. 
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V. THE MAXIMUM LIKELIHOOD ESTIMATE WITH FAILURE 
DISCOUNTING IN THE TEST-FIND-TEST SCENARIO 



A. DESCRIPTION 

The conventional maximum likelihood estimate of 
reliability is: 

R = # of successes/# of trials (5.1) 

In order to use the estimator however, "...constant 
reliability, R, is required for each trial. Because the 
reliability at each phase Rj^ may not be constant with k; 
only the test data from the phase of interest may be used to 
estimate reliability" [Ref. l:p. 22]. An excellent detailed 
description of failure discounting as applied to the test- 
fix-test scenario is contained in Ref . l:pp. 13-21. This 
section primarily addresses the process when employed in 
conjunction with the test-find-test scenario. 

The straight percent discounting method attempts to 
reflect improved system reliability by removing a fraction, 
F, of a failure's weight at an interval of every N trials. 
This sequence of N successful trials must occur in a follow- 
on phase of testing, after the components which were "found" 
to be causes of failure are modified or repaired. The idea 
here is that accumulation of subsequent testing without 
repeat of a failure for the same cause increases the 
confidence that the failure cause has been removed; 
therefore its weight in further estimations should be 
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diminished. When using this method, the data is first 
adjusted and the model is applied to the modified data. The 
success of the model is wholly dependent on the correct 
adjustment of the data and hence the parameters chosen in 
calculating these adjustments. The equation used to compute 
the current value of an adjusted failure is: 

Adjusted failure int (M/N) ( 53 ) 

where M is the number of successful sequential trials in 
follow-on phases for failure cause j and F and N are as 
above. 

B. EXAMPLE 

As an illustration of the method, consider the results 
of the following example. Suppose, as in our previous 

example, testing in each phase continues until the 
occurrence of the third failure. Further suppose the 
results of two phases of testing are as in Table 5. 

TABLE 5. MLEWD EXAMPLE DATA 



Phase 


# Failures 


# Trials to 3rd Failure 


# Successes 


1 


3 


5 


5-3 = 2 


2 


3 


6 


6-3 = 3 
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Then applying equation (5.1), the maximum likelihood 
estimates without failure discounting are: 

= 2/5 = .4 ; R2 = 3/6 = .5 (5.3) 

As one can see, the information obtained from the first 
five trials is disregarded in the calculation of the 
estimate of the reliability in the second phase. Continuing 
in this manner, the data in the first eleven trials will be 
ignored in computation of the third phase reliability. This 
practice is extremely inefficient by any standards. 

In the test-fix-test scenario, [Ref.l:pp. 15-16] it is 
assumed that a design "fix" is implemented after each 
failure; in a test-find-test scenario, no modifications are 
assumed until a given number of failures is observed. The 
causes of the failures are merely noted for subsequent 
action. Since there is no justification for applying the 
discounting method until follow-on phases prove a cause 
corrected, different values of weights are assigned to the 
same sequence of successes and failures in these two 
scenarios. Tables 6 and 7 illustrate the application of the 
discounting method to both situations. In the calculations, 
the parameters N=3 and F=.25 are arbitrarily chosen. 
Unfortunately, in practice, they are often likewise 
selected. 
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TABLE 6 



TEST-FIX-TEST DATA WITH DISCOUNTING 



FAILURE CAUSE 




Notation: For each failure cause the attribute 

S = Success F = Failure 

is listed along with the corresponding weight of that 
attribute • 



TABLE 7 



TEST-FIND-TEST DATA WITH DISCOUNTING 

FAILURE CAUSE 
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Testing in the Test-Find-Test scenario terminates after 
the 3rd failure. The straight percent discounting 
parameters are: 

N = Discount Interval = 3 

F = Fraction the failure is reduced = .25 

Using the test-fix-test discounted data and MLEWD model 
produces the results shown in Table 8 upon completion of the 
eleventh trial . 

TABLE 8. RESULTS OF THE MLEWD MODEL APPLIED TO THE 
TEST-FIX-TEST DATA 



PHASE = 
FAILURE # 


CAUSE 


ADJ. FAILURE 


ADJ. TRIALS 


1 


D 


.75 


1 .T5 = 1..1533 




B 


1 


1 1.0 = 1.0000 


3 


A 


.5625 


3 .5625 = 5.3333 


4 


C 


.75 


2 .75 = 2.6667 


5 


D 


.75 


1 .75 = 1.3333 


6 


B 


1 


II 

q 

r*-, 



TADJ. TRIALS = 14.6665. 

Adj. trials = (the observed = of trials since the previous failure)/ 
( adjusted failure value). 



and 



^ iY ADJUSTED TRIALS) - {TOTAL i^OF FAILURES) 

(> ADJUSTED TRIALS) 



.5909 



(5.4) 
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AS one can see, this estimate is higher than that 
predicted after the eleventh trial (second phase; eg. 5.3) 
using the MLE without discounting. Since it takes into 
account 5 additional data points it is considered to be a 
superior estimate ... provided the parameters were correctly 

chosen. 

The results of the test-f ind-test-scenario after the 
eleventh trial are as shown in Table 9. 

TABLE 9. RESULTS OF THE MLEWD MODEL APPLIED TO THE 
TEST-FIND-TEST DATA 



PHASE 


CAUSE 


ADJ. FAILURE 


AD.I. TRIALS 


I 


D 


1.0 


1/1.0 = 1.0000 


1 


B 


1 


I'l.O = 1.0000 


1 


A 


.5625 


3 .5625 = 5.3333 


2 


C 


1.0 


2T.0 = 2.0000 


2 


D 


1.0 


1 1.0 = 1.0000 


2 


B 


1 


3 1.0 = 3.00U0 



YADJ. TRIALS = 13.3333. 



4 = .5499 

As a comparison, the reliability after the second phase 
of testing in the test-fix-test scenario would be zero since 

there were two failures in two trials. 

The results seem to'indicate that the estimate of 
reliability of a system or equipment will grow more quickly 
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appeal as it seems likely that if failure causes are 
corrected immediately they are else like to be responsible 
for further failure. These results additionally demonstrate 
that the use of discounting has the desirable characteristic 
of producing higher estimates of reliability than the 
conservatively biased MLE without discounting [Ref. 3:p. 
34] . 
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VI. MODEL COMPARISONS 



A. INTRODUCTION 

In sections B through F of this chapter, the performance 
of the unweighted exponential regression model will be 
compared to that of the weighted model used in conjunction 
with weighting methods one and two. The employment of 
method three will be discussed in Section G. The 
performance of the MLEWD will additionally be discussed 
where appropriate. For the purposes of this evaluation 
patterns VI, VII and VIII will be discussed as a group in 
the category. Constant Reliability. Patterns IV and V will 
be similarly grouped under the title Rapid Reliability 
Growth, while pattern I, Convex Reliability, pattern II, 
Decreasing Reliability and pattern III, Intermittent 
Reliability will be evaluated individually. All discussion 
is in regard to test-fix-test scenarios to facilitate 
reference to previous work. 

B. CONSTANT RELIABILITY PATTERNS 

Figure 10 [Ref. 2:p. 32] illustrates the performance of 
the MLEWD model against representative constant reliability 
pattern VI. No failure discounting was implemented in this 
case. The MLEWD model performed better than either the 
weighted or unweighted regression models in all three 
constant reliability patterns, although its superiority 
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diminished with lower levels of constant reliability. This 
was expected since, as indicated in Chapter V, the 
underlying assumption of the MLE is one of constant 
reliability. 

Figure 5 is an illustration of the unweighted and 
weighted exponential regression models also against constant 
reliability pattern VI. Regardless of the level of 
constancy, the unweighted model consistently underestimated 
reliability while the weighted model overestimated it. In 
all patterns, the model weighted using method 2 (i.e., more 
weight to the most recent phase) produced better results 
than when method 1 (i.e., more weight to less variable 
phases) was employed. 

C. RAPID RELIABILITY GROWTH 

Figure 6 depicts the performance of the models when 
applied to Pattern V. This pattern is one in which "actual" 
reliability achieves its greatest value (approximately .9) 
early in the testing process (phase 3) and remains there 
throughout the remainder of the testing. Pattern IV is, in 
all ways, similar to pattern V with the difference being 
that it achieves the value of approximately .99 in phase 3. 
Both versions of the weighted model and the unweighted model 
underestimated reliability during its growth phases although 
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CONSTANT RELIABILITY PATTERN 




Figure 5. Weighting Methods One and Two and the Constant 
Reliability Pattern 

RAPID RELIABILITY GROWTH PATTERN 




Figure 6. Weighting Methods One and Two and the Rapid 
Reliability Growth Pattern 
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the weighted models more accurately tracked the process. As 
with constant reliability, the model weighted by Method 2 
outperformed the model when weighting Method 1 was employed. 
Also similar to constant reliability performance results, 
both versions of the weighted model produced higher 
estimates than the unweighted model over all phases. Given 
this underlying pattern, if one is concerned with tracking 
the initial phases of testing, the weighted model in 
conjunction with Method 2 is recommended as it best captures 
the growth portion of the development process. The 
unweighted model in this case seems to best capture the 
constant reliability portion of the development although its 
superiority is minimal. Figures 13 and 14 [Ref. 2:pp. 36- 
37] illustrate the performance of the MLE both with and 
without discounting. In both cases, the exponential 
regression model outperforms the MLE over all phases. 

D. CONVEX RELIABILITY GROWTH 

This pattern illustrated the most dramatic support for 
application of weighting to the exponential regression 
model. Figure 7 depicts the substantial improvement in 
reliability estimation when weighted Method 2 is employed. 
Use of Method 1 also produced superior results to the 
unweighted model although these were not as dramatic, 
particularly in the latter phases of testing. Both versions 
of the weighted model seemed to capture the pattern of 
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CONVEX RELIABILITY GROWTH PATTERN 




Figure 7. Weighting Methods One and Two and the Convex 
Reliability Growth Pattern 

growth as early as phase 2. No earlier improvements can be 
anticipated as all three versions of the regression model 
require at least two phase of testing to perform a linear 
regression. 

Figure 23 [Ref. 2:p. 46] demonstrates the poor 

performance of the MLE when no discounting is applied. 



39 



These results are not surprising as the underlying 
assumption of constant reliability in the case of the MLE is 
violated here. A dramatic improvement as seen in Figure 24 
[Ref. 2:p. 47] when discount parameters, F=.75 and N=6 are 
applied to the data. The parameters were selected after a 
good deal of trial and error simulation, however, and no 
explanation as to why this choice of parameters was 
effective can currently be offered. Since the luxury of 
simulation is not often available in practice, and the 
weighted model in conjunction with Method 2 outperforms the 
other models, it is recommended for use while the underlying 
reliability pattern is suspected to be convexly increasing. 
E. DECREASING RELIABILITY 

Figure 8 indicates that the same general previous 
observations can be made with regard to the performance of 
the models in the case of decreasing reliability. The 
unweighted regression model produced estimates which are 
everywhere below those of the weighted model regardless of 
the weighting method employed. Additionally, the weighted 
models seem to more accurately capture the changing growth 
pattern. The unweighted model is most accurate when actual 
reliability decreases to a low value however this is 
primarily due to the fact that the unweighted model 
consistently tends to underestimate actual reliability and 
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DECREASING RELIABILITY GROWTH PATTERN 




Figure 8. Weighting Methods One and Two and the 
Decreasing Reliability Growth Pattern 

is less responsive to changes in the actual growth pattern. 
Therefore when actual reliability dipped and the already low 
unweighted model did not respond, the result was a fairly 
accurate estimate. As in previous patterns, the model 
weighted by giving most weight to the most recent phase 
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estimate produces the most responsive results. For this 
reason, this model would be recommended in this case. 

Figures 17 and 18 [Ref. 2 : pp . 40-41] demonstrate the 
inability of the MLE model to capture a decreasing 
reliability pattern whether failure discounting is employed 
or not. 

F. INTERMITTENT RELIABILITY GROWTH 

The Intermittent Reliability Growth pattern (Figure 9) 
can be described as one in which reliability grows fairly 
rapidly, remains constant for several phases and then 
resumes growth. Such a pattern would occur if intended 
improvements to a developing system were ineffective for a 
period of time before discovery of beneficial change was 
implemented. As in all previous patterns, the unweighted 
model produced more conservative estimates than both 
versions of the weighted model. Also as before, the model 
weighted by Method 2 was most responsive to changes in the 
reliability pattern. This is particularly evident in phases 
5 and 6 where it dips below the Method 1 estimate in an 
attempt to recognize the period of constant reliability. 
All three versions of the exponential regression model 
outperformed the MLE model regardless of whether failure 
discounting was employed. Figures 20 and 21 [Ref. 2 :pp. 43- 
44] reflect this. 
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INTERMITTENT RELIABILITY GROWTH PATTERN 




Figure 9. Weighting Methods One and Two and the 
Intermittent Reliability Growth Pattern 



G. THE WEIGHTED MODEL IN CONJUNCTION WITH METHOD 3 
WEIGHTING 

Weighting Method 3 is the application of user selected 
weights to the exponential regression model. Figures 10 - 
14 corresponding to representative examples of constant, 
rapidly increasing, convex, decreasing, and intermittent 
reliability respectively, demonstrate the sensitivity of the 
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model to the choice of weights. In every case the weighted 
regression model identically tracks the unweighted model 
until the user influences the process by applying a weight 
to a phase. For example, when the case 1 weights are 
employed, the user has given equal weight to the first seven 
phases of testing. This is equivalent to not weighting the 
model at all. Therefore, one would expect the unweighted 
model to produce identical results to this version of the 
weighted model. That this is in fact occurring, is most 
clearly seen in Figure 10. Tracking the curve labeled "1", 
(corresponding to case 1, Weights Table 3), one can see it 
is identical to the unweighted growth curve until Phase 7 
after which it abruptly increases. Similar departures from 
the unweighted curve are seen for the other four cases of 
user selected weights. These results are welcome in that 
they substantiate the programming of the weighting process. 
They are disappointing however, in that, as with the MLEWD 
model, the response of the exponential regression model 
seems to be very sensitive to the choice of weights. 
Unlike the MLEWD model though, there are two systematic 
means of selecting weights which consistently produce 
accurate results (i.e., methods 1 and 2). It is recommended 
that these two weighting methods be employed unless the user 
is intimately familiar with the development process and has 
the unique ability to reflect his knowledge in the choice of 
weights . 
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CONSTANT RELIABILITY GROWTH PATTERN 



o 




Figure 10. Method Three Weights and the Constant 
Reliability Growth Pattern 




Figure 11. Method Three Weights and the Rapid 
Reliability Growth Pattern 
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CONVEX RELIABILITY GROWTH PATTERN 




Figure 12. Method Three Weights and the Convex 
Reliability Growth Pattern 




Figure 13. Method Three Weights and the Decreasing 
Reliability Growth Pattern 
46 



RELIABILITY 

0.6 



INTERMITTENT RELIABILITY GROWTH PATTERN 




Fiqure 14. Method Three Weights and the 

Intermittent Reliability Growth 
Pattern 
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VII. SUMMARY, CONCLUSIONS, RECOMMENDATIONS 



A . SUMMARY 

Captain Drake created a FORTRAN program which simulated 
actual reliability patterns and evaluated the performance of 
three reliability growth models against these patterns. His 
work indicated that the MLEWD model was effective against 
all patterns when used in conjunction with the correct 
discounting parameters. He also found that the exponential 
regression model likewise tracked all actual patterns and 
was relatively insensitive to the discounting parameters. 
This model, was however, more variable than the MLEWD model. 

Captain Chandler modified Captain Drake's program to 
produce several additional actual reliability patterns and 
subsequently evaluated both models against these patterns. 
His results demonstrated that the MLEWD model was incapable 
of tracking declining reliability and was only superior to 
the exponential regression model when actual reliability was 
constant; a phenomenon which, by design, rarely occurs in 
the development process. 

The objective of this thesis was to improve the 
exponential regression model by applying weights to the 
model at each phase. The intent was to produce an accurate 
means of estimating changing reliability without burdening 
the user with selecting discounting parameters. 
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Three methods of weighting were derived. The first 
method gives low weight to the most variable estimates of 
reliability, while the second weights most recent results 
more heavily. With both these methods, the FORTRAN program 
systematically computes these weights at each phase. In 
method 3, the user is required to input values of weights 
for each phase of testing. Use of this method produced 
results which were very sensitive to the weights chosen and 
presented the user with the dilemma of selecting critical 
input parameters without any guidance. Since this is the 
very situation which the employment of weights was derived 
to avoid, this method is not recommended. 

B. CONCLUSIONS 

Regardless of the underlying actual reliability pattern 
or the method of weighting, the weighted exponential 
regression model produced higher estimates of reliability 
than the unweighted model. Since the unweighted model 
generally underestimates actual reliability this is 
considered to be a desirable characteristic. Both the 
weighted and unweighted model produced results comparable or 
superior to the MLEWD model against all variations of the 
actual reliability pattern with the exception of the case of 
constant reliability. In this case, the MLEWD would be 
recommended . 
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In all other cases, the regression model weighted by 
Method 2 seems to be most responsive to changes in actual 
reliability. Additionally, it produced more accurate 
estimates than the model weighted by either Methods 1 or 3. 
In the case of convex reliability Method 2 weights are a 
clear choice over the unweighted model as well. When a 
period of decreasing reliability is suspected in the 
development process, the Method 2 weights are also 
recommended unless one is concerned with the magnitude of 
the decrease as opposed to the trend of the pattern. In 
this case the more conservative unweighted model is 
recommended. When estimating the reliability of the 
intermittent and rapidly increasing patterns, the Method 2 
weighted model most accurately captured the growth phases of 
reliability but were less accurate than the unweighted model 
in the final phase of testing although the differences in 
estimates here were minimal. Based on these conclusions, it 
is recommended that if a user has little or no knowledge of 
the actual reliability growth of a developing system or 
equipment, he would be wise to select the exponential 
regression model in conjunction with Method 2 weights to 
estimate the reliability growth curve. 

As a final conclusion to this study, implementing the 

ic 

unbiased estimator with minimum variance, Yj^ in all versions 
of the exponential regression model resulted in an improved 
reliability estimate over the models in conjunction with Yj^. 
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C. RECOMMENDATIONS FOR FURTHER STUDY 



The following are recommendations for further study. 

1. Since a large segment of the commercial world collects 
its reliability data in terms of mean time between 
failures (MTBF's) it is highly recommended that a 
simulation to analyze continuous data be derived. 

2. Presently, straight percent failure discounting has 
shown some promising results when the discounting 
parameters are correctly chosen. Development of 
guidelines for selecting these parameters would 
certainly be of value in the employment of the MLEWD 
model. A difficulty in studying these parameters is 
that they cannot be altered from phase to phase or 
varied with failure cause in the current version of 
the program . 

3. As with failure discounting, the employment of user 
selected weights might prove beneficial if guidelines 
could be developed to aid in their selection. A more 
extensive analysis of this weighting method is 
recommended . 
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APPENDIX A 



USER'S GUIDE TO DISCRETE RELIABILITY GROWTH ( DRG ) 

1. Introduction 

2. The DRG Exec File 

3. The Input Data File 
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APPENDIX A: USER'S GUIDE TO DISCRETE 

RELIABILITY GROWTH (DRG) 

1 . Introduction 

In order to use the Fortran program, DRG, the user must 
possess three files: 

1. Input_ data A1 

2. DRG Fortran Al 

3. DRG Exec Al 

A sample of each of these files is contained in Appendix 
B along with sample output. The input file and the exec 
file may be tailored to the user's needs. In its current 
form, the exec file produces a large degree of intermediate 
calculations for both the DRG Fortran program and Captain 
Chandler's program, JIMC Fortran Al . A copy of this program 
is also contained in Appendix B. These calculations may not 
be of interest to the user and may be eliminated with no 
detrimental effect to the program. A detailed explanation 
of each file definition is contained in section 2 of this 
appendix. 

The difference between the two versions of the Program 

is that JIMC Fortrans retains the ability to employ 

— ★ 

the weighted average estimate model and uses Y vice Y in 
the calculations of all versions of the exponential 
regression model. (see chapter III for a discussion of the 
motivation for altering this parameter which is discussed in 
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Section 3 para 12. In all other ways, the input files for 
these two programs are identical . 

2. The DRG Exec File 

This file contains all the requisite file definitions 
and commands to run either DRG Fortran A1 or JIMC Fortran 
Al. To indicate which of these programs is desired, the 
user must alter the third line of the file so that it reads, 
&FN=DRG or &FN=JIM respectively. Once this has been done, 
the user simply types DRG or JIMC while in CMS to execute 
the program. 

Of the remaining lines, only those beginning with the 
word FILEDEF should be altered or eliminated. All other 
lines contain commands pertaining to the execution of the 
program. The first FILEDEF line currently reads FILEDEF 10 
DISK INPUTI DATA Al . The "1" in the filename indicates that 
this data file produces actual reliability growth pattern 1 
described in this thesis. Eight patterns and input files 
have been established for each program. These input files 
may be duplicated as described in Section 3 using the 
corresponding pattern parameters in Ref. 2 if the user does 
not have them readily available. This FILEDEF may not be 
eliminated from the exec file, however the "1" may be 
changed to any integer up to and including "8" for use with 
DRG Fortran Al. The integers 9 through 16 correspond to the 
same patterns when executing JIMC Fortran Al . This 
numbering convention is adhered to with the remaining 
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FILEDEFs. Since a filename is only permitted to be six 
characters in length, whenever necessary, it is truncated to 
include the integers discussed. For example, input 10 is an 
unacceptable filename therefore the "T" is dropped from the 
word "input" to allow "inpul0." 

The next two FILEDEFs have the filenames Al and A9 
respectively. The first corresponds to the "actual" 
reliability pattern (pattern 1) produced by the DRG program 
while the second corresponds to the same pattern when JIMC 
Fortran Al is run. These files should contain identical 
output. Both are kept for convenience so that if for 
instance, the user runs DRG with pattern 1 input and then 
JIMC with pattern 2 input he will still retain the initial 
"actual" reliability file, i.e., it will not be overwritten. 
These two files were both devised solely for plotting 
purposes and may be eliminated with no disruption to the 
program. The information contained in these files is 
duplicated in PRELIAB and JRELIAB Listing respectively. 

PRELIAB LISTING Al contains the primary comprehensive 
output for the DRG Fortran program while JRELIAB contains 
similar information for JIMC Fortran Al . These files 
consist of a model parameter summary, a comparison of the 
mean predicted reliability for each model to the actual, the 
estimate standard deviation and a 95% confidence interval 
for each model. Also contained in these files is a 
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recapitulation of results so that each model may be compared 
to all others. These FILEDEFs should not be eliminated from 
the exec file. 

JThesis and PThesis are the filenames corresponding to 
the intermediate summary output produced by JIMC and DRG 
Fortran respectively. These files contain a phase by phase 
listing of failures and failure causes for up to 5 
simulations. This output may be "turned on (or off)" by 
placing a 1 (or 0)" in the appropriate line of the input 
file. See section 3 para 12 for details. These FILEDEFs 
should not be eliminated from the exec file. 

PMATRIXA and JMATRIXA contain the working "A" matrix 
for each of the programs. The working A matrix contains 
such information as the probability of success for each 
failure cause, the number of trials to failure for each 
cause, the system cause of failure, the phase number, and, 
in the case of JMATRIXA, the adjusted number of trials and 
adjusted number of failures. These FILEDEFs may be 
eliminated if desired. 

PREGMAT and JREGMAT contain the parameters computed in 
the REG, WREGl, WREG2, and WREG3 arrays as described on page 
2 of the program heading, (Appendix B, section 3). These 
mat ices pertain to the exponential reg ression and Weighted 
exponential regression models. These files are not 
required. 
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Ystar LISTING is produced only by DRG Fortran A1 . It 
currently lists (arbitarily) the first 219 values for the 
tenth phase of testing. It is included merely as a 
verification of the implementation of this parameter. This 
FILEDEF may be eliminated. 

TRIALS DATA A1 is also only produced by DRG Fortran A1 . 
It lists the adjusted trials to failure for each phase. In 
addition, the weights computed for use by methods 1 and 2 
are contained here. This FILEDEF is not required. 

EST OUT A1 is produced only by JIMC Fortran A1 and 
contains each phase estimate generated by the Woods Weighted 
Average Estimate Model for up to the first 5 simulations 
run. The remaining FILEDEFs with File Type "OUT" contain 
similar information generated by the filename models. A 
file name beginning with the letter "J" corresponds to 
output produced by JIMC Fortran Al, while a file name 
starting with the letter "P" originated from DRG Fortran Al. 
These FILEDEFs should be retained. 

All of the remaining FILEDEFs correspond to files 
generated solely for plotting purposes. Each contains the 
mean reliability estimates of each phase in the Filetype 
"NUM" file and the corresponding standard deviation of these 
estimates in the Filetype "SDV" file. The information in 
these Files is duplicated in PRELIAB and JRELIAB LISTING. If 
the user is not interested in graphing results, these 24 
FILEDEFs may be eliminated. As explained earlier, the 



57 



integers in the filenames of these files correspond to the 
underlying actual reliability growth pattern. 

3. The Input Data File 

These instructions should be used in conjunction with a 
sample input file (See Appendix B, section 2) to aid in 
understanding. The easiest method for preparing an input 
file is to edit an existing input file. Due to the 
formatting of the "READ" statements in the simulation, it is 
imperative that all inputs be entered in the correct 
sequence. The simulation is built to read input from device 
number 10. 

The following steps should be allowed in order to 
produce an acceptable input file. All entries must be on a 
separate line in the input file although they may be 
anywhere in the line as long as the data entry is the first 
item encountered. 



1. Determine how many failure causes will be allowed in 
this simulation. This number must be an integer 
greater than or equal to one. There is no set limit 
on how large a number is possible. The capacity of 
the machine on which the program is being run will 
have some effect. Enter this number on the top line. 

2. Determine how many test phases are desired. Again, 
this number must be an integer greater than one. 
Enter this number on line 2. 

3. If the fixed phase reliability mode is desired than 
enter a 1 on line 3. If the constrained random growth 
mode is required than enter a 0 on line 3. 
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4 . 



Enter the number of failures that will be allowed in 
Phase 1 on line 4. This number must also be an 
integer. A test-fix-test scenario would have one 
failure per phase, for example. Repeat this process 
on successive lines until you have entered the number 
of failures allowed for each phase that will be 
tested. Remember, if you defined the test as being 10 
phases, then you should have 10 separate entries, each 
on its own line, for this step. 

5. If the constrained random growth option is selected 
(you should have entered a 0 at line 3) then you must 
now enter the probability of success due to each 
failure cause for the first phase. This step should 
contain as many lines as the number of failure causes 
you identified in Step 1. The probability of success 
due to a failure cause is simply the probability of 
the failure cause occurring subtracted from one. 

6. If the fixed phase reliability option is selected (you 
should have entered a 1 at line 3) then the process of 
entering the probabilities of success is a little 
more involved. Fixed phase reliability means that one 
is controlling the actual system reliability at each 
phase of the test. Since this system reliability is 
merely the product of the probabilities of success of 
all the failure causes at each phase then the user 
must enter these probabilities. These probabilities 
are entered by failure cause for each phase. If there 
are two failure causes and 3 phases in a test then the 
first entry will be the probability of success due to 
failure cause 1 in phase 1; the next entry will be the 
probability of success due to failure cause 1 in phase 
2 and then failure cause 1 in phase 3. After failure 
cause one has been entered for all three phase then 
failure cause 2's probabilities of success should be 
entered by by phase. Remember that each data entry 
must be on its own separate line. If X is the number 
of failure causes and Y is the number of test phase 
then this step should result in XY total data input 
lines. This number does not alter if more than one 
failure is permitted in a phase. 
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7 . 



The next item that must be entered concerns use of the 
standard failure discount method. Even if you desire 
to use the C.L. discount method or if you do not 
desire to discount previous failures at all you must 
still place a value in this line. This line requires 
you to enter the discount interval or the number of 
successful trials that must occur between applications 
of the standard discount method. This number should 
be an integer. 

8. The next item required is the discount fraction. This 
value also applies to the standard discount method. 
If you do not want to discount previous failures then 
your should enter 0.0 on this line. If you desire to 
discount previous failures by 50 percent each time the 
discount method is applied then you should enter 0.50 
on this line. If you are using the C.L. discount 
method, you must still enter a value here, although it 
will not be used by the models in estimating 
reliability. 

9 . The next required entry is the random number seed for 
random number generator. Any number greater than zero 
and less than 2^^ will suffice. 

10. The FRIMP must be entered next. This value is only 
used if the constrained random growth mode is selected 
at line 3. This number must be between 0 and 1 and 
represents the fraction by which reliability will 
improve from phase to phase. It is applied to the 
probabilities of success of the failure causes each 
time they cause system failure in a phase. This 
method is intended to represent repairs or 
improvements in the system during the test. An entry 
of 0 will result in constant reliability while entries 
close to one will result in rapid reliability growth. 
Even if you are using the fixed phase reliability 
option you must still enter a value here although it 
will not play any role in the simulation. 

11. The next item that must be entered is the number of 
replications desired. For the purposes of the thesis 
associated with this paper 500 replications were done 
for each reliability growth pattern. 
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12. The next 7 (8, for use with JIMC Fortran A1 ) lines of 
required input deal with intermediate and specified 
output. They are all binary (0 or 1) option 
statements. The first addresses intermediate output. 
If intermediate phase by phase data is desired then 
the user should enter a 1 here. If intermediate 
output is not desired then enter a 0. This output 
will get voluminous rather quickly, particularly if a 
large number of relications is requested. The 
remaining options deal with the estimates generated 
by each model at each phase. The final output (which 
will be given regardless of the options chosen here) 
represents the average of all these estimates. If 
phase estimates are desired then a 1 should be 
entered; if not then enter a 0. The order of entry is 
weighted average model (for JIMC Fortran program 
only), MLE with failure discounting model, MLE single 
phase model, exponential regression model, and 
weighted exponential regression methods 1 through 3. 
Again, each data entry must be on a separate line. 

13. The next required entry is the failure discounting 
option. If you desire to use the standard discount 
method then enter a 1. If you desire not to discount 
at all enter a 1 and make sure you have entered 0.0 as 
instructed at Step 8 above. If you desire to use the 
C.L. method of failure discounting then enter a 2 in 
this space. 

14. Next enter the value of the C.L. method parameter (the 
confidence interval ) . This number should be between 0 
and 1 and will typically be in the range of 0.8 to 
0.99. 

15. The next item that must be provided is the C.L. 
discount interval. Normally, the C.L. method does not 
use a discount interval (defined similar to the 
definition of a discount interval for the standard 
discount method) but in the majority of cases applying 
this method with some type of specified interval will 
lead to better results. These last two items must be 
entered even if one is using a alternate method of 
discounting. The other method will be applied but 
values are required due to the particular formatting 
of the simulation. 

16. The final entries that must be made are the weights to 
be accorded each phase of development. Enter any 
fraction between zero and one for each phase on a 
separate line. Remember the sum of these fractions 
must equal one. These entries are required regardless 
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of whether the user desires to apply weighting Method 
3 or not. Note: Giving equal weight to all phase is 
equivalent to not weighting and will produce results 
identical to the exponential regression model . 



If the input file is established consistent with the 
format outlined above and with the sample input file 
provided in Section 2 of Appendix B then there should be no 
problem in obtaining results. 
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APPENDIX B: FILES AND PROGRAMS 



1. SAMPLE EXEC FILE 

2. SAMPLE INPUT DATA FILE 

3. DRG FORTRAN PROGRAM 

4. JIMC FORTRAN PROGRAM 

5. SAMPLE OUTPUT FROM DRG FORTRAN A1 
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1. SAMPLE EXEC FILE 



&TRACE OFF 



**************************************************************************** 



&FN = DRG 
&FN1 = RELIAB 

&TYPE Do you need to compile your program ? (Y) 

&READ VAR &R_COMPILE 

&IF &R_COMPILE NE Y &G0T0 -RUN 

-H FORTVS &FN 

&IF &RC EQ 0 &G0T0 -RUN 

&TYPE Your program did not compile; check for errors. 
&TYPE Do you wish to view the program LISTING file? (Y) 
&READ VAR &RSP1 

&IF &RSP1 EQ Y BROWSE &FN LISTING A 

&TYPE Do you wish to XEDIT the program file? (Y) 

&READ VAR &RESP1 

&IF &RESP1 NE Y &EXIT 1 

&COMMAND XEDIT &FN FORTRAN A 

&TYPE Do you wish to run the program again? (Y) 

&READ VAR &RESP2 

&IF &RESP2 EQ Y &G0T0 -H 

&EXIT 1 

-RUN 

FILEDEF 10 DISK INPUTl DATA A1 

FILEDEF 82 DISK A1 NUM A1 

FILEDEF 84 DISK A9 NUM A1 

FILEDEF 30 DISK JRELIAB LISTING A1 (LRECL, 133 

FILEDEF 35 DISK PRELIAB LISTING Al (LRECL 133 

FILEDEF 20 DISK JTHESIS OUT Al 

FILEDEF 25 DISK PTHESIS OUT Al 

FILEDEF 81 DISK PMATRIXA LISTING (LRECL 133 

FILEDEF 83 DISK JMATRIXA LISTING (LRECL 133 

FILEDEF 87 DISK PREGHAT DATA Al 

FILEDEF 88 DISK JREGMAT DATA Al 

FILEDEF 90 DISK YSTAR LISTING (LRECL 133 

FILEDEF 89 DISK TRIALS DATA Al 

FILEDEF 40 DISK EST OUT Al 

FILEDEF 50 DISK JMLEWD OUT Al 

FILEDEF 55 DISK PMLEWD OUT Al 

FILEDEF 60 DISK JMLESP OUT Al 

FILEDEF 65 DISK PMLESP OUT Al 

FILEDEF 70 DISK JREGEST OUT Al 

FILEDEF 75 DISK PREGEST OUT Al 

FILEDEF 15 DISK PWRESl OUT Al 

FILEDEF 39 DISK PWRES2 OUT Al 

FILEDEF 49 DISK PWRES3 OUT Al 

FILEDEF 16 DISK JVRESl OUT Al 

FILEDEF 38 DISK JVRES2 OUT Al 

FILEDEF 48 DISK JW'RES3 OUT Al 

FILEDEF 52 DISK MLEWDl NUM Al 

FILEDEF 51 DISK MLEWDl SDV Al 

FILEDEF 54 DISK MLEWD9 NUM Al 
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FILEDEF 


53 


DISK 


MLEWD9 SDV A1 


FILEDEF 


72 


DISK 


REG8 


NUM 


A1 


FILEDEF 


71 


DISK 


REGS 


SDV 


A1 


FILEDEF 


74 


DISK 


REG 16 NUM A1 


FILEDEF 


73 


DISK 


REG16 SDV A1 


FILEDEF 


77 


DISK 


MlPl 


NUM 


A1 


FILEDEF 


76 


DISK 


MlPl 


SDV 


A1 


FILEDEF 


79 


DISK 


M2P1 


NUM 


A1 


FILEDEF 


78 


DISK 


M2P1 


SDV 


A1 


FILEDEF 


92 


DISK 


M3P1 


NUM 


A1 


FILEDEF 


91 


DISK 


M3P1 


SDV 


A1 


FILEDEF 


18 


DISK 


M1P9 


NUM 


A1 


FILEDEF 


17 


DISK 


M1P9 


SDV 


A1 


FILEDEF 


94 


DISK 


M2P9 


NUM 


A1 


FILEDEF 


93 


DISK 


M2P9 


SDV 


A1 


FILEDEF 


96 


DISK 


M3P9 


NUM 


A1 


FILEDEF 


95 


DISK 


M3P9 


SDV 


A1 


FILEDEF 


62 


DISK 


MLESPl NUM A1 


FILEDEF 


61 


DISK 


MLESPl SDV A1 


FILEDEF 


64 


DISK 


MLESP9 NUM A1 


FILEDEF 


63 


DISK 


MLESP9 SDV A1 


FILEDEF 


06 


TERMINAL 







LOAD &FN (START 
&IF ARC EQ 0 ASKIP 9 

ATYPE Your program did not run correctly; check for errors. 
ATYPE Do you wish to XEDIT the program file? (Y) 

AREAD VAR ARESP3 

AIF ARESP3 NE Y AEXIT 2 

ACOMMAN'D XEDIT AFN FORTRAN A 

ATYPE Do you wish to run the program again? (Y) 

AREAD VAR ARESP4 

AIF ARESP4 EQ Y AGOTO -H 

AEXIT 2 

ATYPE YOUR OUTPUT IS IN THE FILE AFNl LISTING A 
ATYPE Do you wish to BROWSE your output? (Y) 

AREAD VAR ARESP 

AIF ARESP EQ Y ACOMMAND BROWSE AFNl LISTING A 
ATYPE Print your output file? (Y) 

AREAD VAR ARESP7 

AIF ARESP7 EQ Y ACOMMAND PRINT AFN LISTING A 
-REDO 

ATYPE Do you wish to XEDIT the program file? (Y/N) 

AREAD VAR ARESP5 

AIF ARESP5 EQ Y XEDIT AFN FORTRAN A 

ATYPE Do you wish to run the program again? (Y) 

AREAD VAR ARESP6 

ARESP56 = ACONCAT OF ARESP5 ARESP6 
AIF ARESP56 EQ YY AGOTO -H 
AIF ARESP6 EQ Y AGOTO -RUN 
AEXIT 
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5 

10 

1 

1 

1 

1 

1 

1 

1 

1 

1 

1 

1 

0. 85 

0 . 86 

0. 90 

0. 91 

0. 93 

0. 95 

0. 97 

0. 99 

0. 99 

0. 998 

0. 84 

0. 85 

0. 87 

0. 90 

0. 92 

0. 95 

0. 97 

0. 99 

0. 99 

0. 998 

0. 83 

0. 84 

0 . 86 

0 . 88 

0. 90 

0. 93 

0. 96 

0. 98 

0. 99 

0. 998 

0. 83 

0. 84 

0. 85 

0. 87 

0. 89 

0. 92 

0. 94 

0. 975 



2. SAMPLE INPUT DATA FILE 



NUMBER OF FAILURE CAUSES 
NUMBER OF PHASES ( NPHASE ) 

FIXED RELIABILITY OPTION ( 1: YES ; 0: NO ) 



NUMBER 


OF 


FAILURES 


IN PHASE 1 










NUMBER 


OF 


FAILURES 


IN PHASE 2 










NUMBER 


OF 


FAILURES 


IN PHASE 3 










NUMBER 


OF 


FAILURES 


IN PHASE 4 










NUMBER 


OF 


FAILURES 


IN PHASE 5 










NUMBER 


OF 


FAILURES 


IN PHASE 6 










NUMBER 


OF 


FAILURES 


IN PHASE 7 










NUMBER 


OF 


FAILURES 


IN PHASE 8 










NUMBER 


OF 


FAILURES 


IN PHASE 9 










NUMBER 


OF 


FAILURES 


IN PHASE 10 








PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


1 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


2 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


3 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


4 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


5 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


6 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


7 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


8 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


9 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


1 


IN 


PHASE 


10 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


1 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


2 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


3 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


4 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


5 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


6 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


7 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


8 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


9 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


2 


IN 


PHASE 


10 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


1 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


2 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


3 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


4 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


5 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


6 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


7 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


8 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


9 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


3 


IN 


PHASE 


10 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


1 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


2 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


3 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


4 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


5 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


6 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


7 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


8 
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0. 99 
0. 998 
0. 81 
0. 83 
0. 84 
0. 86 
0. 89 
0. 91 
0. 94 
0. 961 
0. 99 
0. 998 
1 

0 . 0 

624712. 0 
0. 75 
500 
0 
0 
0 
0 
0 
0 
0 
1 

0. 9 
1 

0. 03 
0. 03 
0. 03 
0. 03 
0. 03 
0. 15 
0. 15 
0. 15 
0 . 2 
0 . 2 



PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


9 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


4 


IN 


PHASE 


10 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


1 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


2 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


3 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


4 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


5 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


6 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


7 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


8 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


9 


PROB. 


OF 


SUCCESS 


FROM 


CAUSE 


5 


IN 


PHASE 


10 


NUMBER 


OF 


TRIALS AFTER 


FAILURE 


BEFORE A 


DISCOUNT IS APPLIED 



FRACTION EACH FAILURE IS DISCOUNTED 

RANDOM NUMBER SEED FOR GGUBFS UNIFORM (0,1) GENERATOR 
FRACTION RELIABILITY IMPROVES AFTER FAILING IN A PHASE 
NUMBER OF DESIRED REPETITIONS FOR THE SIMULATION 
INTERMEDIATE INPUT OPTION(l; INT. OUT; 0: NO INT. OUTPUT) 
SAVE ALL MLE W/ DISCOUNTING ESTIMATES (1: YES; 0: NO ) 

SAVE ALL MLE SINGLE PHASE ESTIMATES (1: YES; 0: NO ) 

SAVE ALL UNVT’D REGRESSION ESTIMATES (1: YES; 0: NO ) 

SAVE ALL METHOD 1 WT’D REG. ESTIMATES (1: YES; 0: NO ) 

SAVE ALL METHOD 2 WT’D REG. ESTIMATES (1: YES; 0: NO ) 

SAVE ALL METHOD 3 WT’D REG. ESTIMATES (1: YES; 0: NO ) 

DISCOUNTING OPTION (1: STRAIGHT % ; 2: LLOYD METHOD) 

PERCENT C. I. FOR C. L. DISCOUNTING METHOD(MUST HAVE A VALUE 
C. L. DISCOUNT INTERVAL 
WEIGHT FOR PHASE 1 
WEIGHT FOR PHASE 2 
WEIGHT FOR PHASE 1 
WEIGHT FOR PHASE 2 
WEIGHT FOR PHASE 1 
WEIGHT FOR PHASE 2 
WEIGHT FOR PHASE 1 
WEIGHT FOR PHASE 8 
WEIGHT FOR PHASE 9 
WEIGHT FOR PHASE 10 
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3. DRG FORTRAN PROGRAM 






* DISCRETE RELIABILITY GROWTH SIMULATION 

it 



it 

it 



* PROGRAMMED BY : JAMES E. DRAKE, 

* JAMES D. CHANDLER, 

* AND PAM A. MARKIEWICZ 

it 



it 



it 



it 



LAST MODIFIED 11 JUN 1988 

THE FOLLOWING EXTERNAL FILES ARE USED BY THE PROGRAM 

INPUT : DATA AND PARAMETER INPUT FILE (DEVICE # 10) 

THESIS : OUTPUT FILE CONTAINING INTERMEDIATE COMPUTATIONS 
(DEVICE # 25) 

RELIAB: OUTPUT FILE CONTAINING FINAL RESULTS OF THE SIMULATION 
(DEVICE 35) 

MLEWD : OUTPUT FILE CONTAINING MLE ESTIMATES USING DISCOUNTING 
FOR EACH PHASE AND EACH REPLICATION 
(DEVICE il: 55) 

MLESP : OUTPUT FILE CONTAINING MLE ESTIMATE FOR EACH SINGLE PHASE 
AND ALL REPLICATIONS USING NO DISCOUNTING 
(DEVICE # 65) 

REGEST : OUTPUT FILE CONTAINING EACH PHASE ESTIMATE FOR EACH 
REPLICATION OF THE EXPONENTIAL REGRESSION ESTIMATE 
(DEVICE 75) 

WRESTl : OUTPUT FILE CONTAINING EACH PHASE ESTIMATE FOR EACH 
REPLICATION OF THE WEIGHTED EXPONENTIAL REGRESSION 
ESTIMATE USING WEIGHTING METHOD 1 (DEVICE ^)15) 

WREST2 : SAME AS ABOVE USING WEIGHTING METHOD 2 (DEVICE j')39) 

WREST3 : SAME AS ABOVE USING WEIGHTING METHOD 3 (DEVICE #49) 



* THE FOLLOWING IS A LIST OF KEY ARRAYS USED IN THE SIMULATION 



y« 

* 

JU 



Vf 

yc 

y« 

y? 

y? 



yc 



yc 

y» 



yc 

Vc 



Vc 

* 

yc 



* A 

* 



* NFAPH 

yc 

* NFCAUS 

yc 

* NTRIAL 
yc 

* TRIALS 

* PHREST 



MAIN WORKING ARRAY CONTAINS PROBABILITY OF SUCCESS FOR 
EACH FAILURE CAUSE, NUMBER OF TRIALS UNTIL FAILURE FOR 
EACH FAILURE CAUSE AND THE SYSTEM, CAUSE OF FAILURE, 
PHASE NUMBER, ADJUSTED NUMBER OF TRIALS AND ADJUSTED 
NUMBER OF FAILURES 

DIMENSION ( ((2*#CAUSES)+6) , (/FAILURES ) 

CONTAINS THE NUMBER OF FAILURES IN EACH PHASE 
DIMENSION (1, (/PHASES) 

BINARY ARRAY USED TO DETERMINE IF A FAILURE OCCURRED IN 
A PHASE 

DIMENSION ( 1, (/FAILURE CAUSES) 

CONTAINS THE NUMBER OF TRIALS SINCE LAST FAILURE OR 
DISCOUNTING FOR EACH FAILURE CAUSE 
DIMENSION ( 1, (/FAILURE CAUSES ) 

CONTAINS THE TOTAL NUMBER OF ADJUSTED TRIALS IN A PHASE 
DIMENSION ( 1, (/ OF PHASES) 

RECORDS THE PHASE ESTIMATE FOR EACH ESTIMATOR WITHIN A 
SINGLE REPLICATION 



* 



* 



* 

* 

* 



* 



* 

yc 

yc 



* 

* 

* 
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* 




DIMENSION (6,f'/PHASES) 


* 


* 




ROW 1 


MLE WITH DISCOUNTING 


* 


* 




ROW 2 


SINGLE PHASE MLE 


k 


* 




ROW 3 


WOODS REGRESSION ESTIMATE (UNWEIGHTED) 


k 


•)V 




ROW 4 


WEIGHTED REGRESSION ESTIMATE (METHOD 1) 


k 


* 




ROW 5 


WEIGHTED REGRESSION ESTIMATE (METHOD 2) 


k 


* 




ROW 6 


WEIGHTED REGRESSION ESTIMATE (METHOD 3) 


k 




AREL 


: CONTAINS ACTUAL SYSTEM RELIABILITY IN EACH PHASE 


k 


Vr 




DIMENSION (1,#PHASES) 


k 


* 


YJK 


: CONTAINS YJK VALUES UP TO 1000 


k 


* 




DIMENSION (1,1000) 


k 


* 


YSTAR 


: CONTAINS THE YK VALUES UP TO 1000 


k 


Vf 




DIMENSION (#PHASES,1000) 


k 


•jV 


VAR 


: CONTAINS THE VARIANCE OF THE ESTIMATE AT EACH PHASE 


k 


•k 




DIMENSION (1, //PHASES) 


k 


k 


SVARl 


: CONTAINS THE RUNNING SUM OF 1/VAR(K) UP TO THE KTH PHASE 


k 


k 




DIMENSION (1, //PHASES) 


k 


k 


SVAR2 


: CONTAINS THE RUNNING SUM OF VAR(K) UP TO THE KTH PHASE 


k 


k 




DIMENSION (1, //PHASES) 


k 


k 


SUMW 


: CONTAINS THE RUNNING SUM OF WEIGHTS UP TO THE KTH PHASE 


k 


k 




DIMENSION (1, //PHASES) 


k 


k 


W 


: CONTAINS THE WEIGHT GIVEN TO EACH RELIABILITY ESTIMATE AT 


k 


k 




EACH PHASE. DIMENSION (3, //PHASES) 


k 


k 

k 




ROW 1: 


WEIGHT = (1/VAR(K))/(SUM OF 1/VAR(K) UP TO 
THE KTH PHASE) 


k 

k 


k 




ROW 2: 


WEIGHT = VAR(K)/(SUM OF VAR(K) UP TO KTH PHASE) 


k 


k 




ROW 3: 


WEIGHT = USER INPUTTED VALUES. 


k 


k 


REG 


: ARRAY USED TO COMPUTE THE EXPONENTIAL REGRESSION ESTIMATE 


k 


k 




DIMENSION (5, //PHASES) 


k 


k 




ROW 1 


K BAR 


k 


k 




ROW 2 


Y BAR 


k 


k 




ROW 3 


Y BAR FOR THE PHASE 


k 


k 




ROW 4 


B HAT 


k 


k 




ROW 5 


A HAT 


k 


k 


WREGl 


; ARRAY USED TO COMPUTE THE WEIGHTED EXPONENTIAL REGRESSION 


k 


k 




ESTIMATE USI.NG METHOD 1 WEIGHTS. DIMENSION (4, //PHASES) 


k 


k 




ROW 1 : 


; X BAR(W) ...(A.K.A. K-BAR(W)) 


k 


k 




ROW 2 : 


: B HAT(W) 




k 




ROW 3 : 


: A HAT(W) 




k 

k 




ROW 4 ; 


: Y BAR(W) 


* 




WREG2 


: SAME AS ABOVE USING METHOD 2 WEIGHTS 




k 


WREG3 


: SAME AS ABOVE USING METHOD 3 WEIGHTS 


k 



* THE REMAINING ARRAYS ARE USED TO COMPUTE THE MEAN AND VARIANCE * 

* OF EACH ESTIMATE AT EACH PHASE. THEY ALL HAVE THE SAME DIMENSIONS * 

* AND STRUCTURE * 

DIMENSION (4,#PHASES) * 

ROW 1 : RUNNING SUM OF ESTIMATES * 

ROW 2 : RUNNING SUM OF SQUARED ESTIMATES * 

ROW 3 : MEAN OF THE ESTIMATES * 

ROW 4 : STANDARD DEVIATION OF THE ESTIMATES * 

Vf 

VALUES FOR THE MLE WITH DISCOUNTING * 

VALUES FOR THE SINGLE PHASE MLE 

VALUES FOR THE EXPONENTIAL REGRESSION ESTIMATE * 

VALUES FOR THE WEIGHTED EXPONENTIAL REGRESSION ESTIMATE * 



* 

* 

* 



MLEWD 

MLESP 

REGEST 

WRESTl 
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* USING METHOD 1 WEIGHTS. * 

* WREST2 : VALUES FOR THE WEIGHTED EXPONENTIAL REGRESSION ESTIMATE * 

* USING METHOD 2 WEIGHTS. * 

* WREST3 : VALUES FOR THE WEIGHTED EXPONENTIAL REGRESSION ESTIMATE * 

* USING METHOD 3 WEIGHTS. 






C DEFINE AND DIMENSION VARIABLES 

PARAMETER (NR=50 ,NC=200) 

INTEGER REP , DI SOPT , FRELOP , LDI , ALD 
REAL*4 MIN 

REAL*8 DSEED , MLE SP , MLE WD , EUL , SUM 1 , SUMS 1 , SUM2 , SUMS 2 , SUM3 
C, SUMS 3, SUM 1 A 

DIMENSION NFAPH(NR) ,A(NR,NC) ,NFCAUS(NR) ,NTRIAL(NR) ,PHREST( 6 ,NR) , 
CMLEWD(4,NR) ,MLESP(4,NR) ,REGEST(4,NR) , AREL(NR) , YJK( 1000) , 
CREG(5,NR) ,YSTAR(50, 1000) ,TADJTP( 1000) ,VAR(NR) ,SVAR1(NR) ,W(3,NR) , 
CWREG1(4,NR) ,WREST1(4,NR) ,WREST2(4,NR) ,WREST3(4,NR) , 

CSVAR2(NR) ,WREG2(4,NR) , WREG3( 4,NR) , SUMW(NR) 

C READ IN THE NUMBER OF CAUSES TO BE USED ( NCAUSE ) AND THE NUMBER 
C OF PHASES ( NPHASE ) IN THE TEST 

READCIO,*) NCAUSE 
READ(10,*) NPHASE 

C CHECK IF FIXED RELIABILITY OPTION IS DESIRED. FIX EULER* S NUMBER. 



READCIO,*) FRELOP 
EUL = 0. 5772156648 

C CREATE VARIABLES FOR THE ROW INDICES OF THE WORKING MATRIX ( A ) 

C IPHASE: PHASE 

C ISYSPR: ACTUAL COMPONENT RELIABILITY 

C INTR: NUMBER OF TRIALS UP TO AND INCLUDING FAILURE 
C IFAILC: CAUSE OF THE FAILURE 

C lADJF: ADJUSTED NUMBER OF FAILURES ED 

C AFTER DISCOUNTING HAS BEEN APPLIED 

C lADJT: ADJUSTED NUMBER OF TRIALS AFTER DISCOUNTING HAS BEEN APPLIED 

IPHASE = (2*NCAUSE)+1 
ISYSPR = IPHASE +1 
INTR = ISYSPR + 1 
IFAILC = INTR + 1 
lADJF = IFAILC + 1 
lADJT = lADJF + 1 



C READ IN THE NUMBER OF FAILURES IN EACH PHASE ( NFAPH(I) ) AND 
C COMPUTE THE TOTAL NUMBER OF FAILURES IN THE TEST ( NFAIL ) 

NFAIL = 0 

DO 10 1=1, NPHASE 
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READCIO,*) NFAPH(I) 

NFAIL = NFAIL + NFAPH(I) 
10 CONTINUE 



C INPUT THE PROBABILITY OF SUCCESS IN A SINGLE TRIAL FOR EACH CAUSE 
C IN EACH PHASE IF FRELOP EQUALS ONE. 

IF (FRELOP .EQ. 1) THEN 
DO 15 I=1,NCAUSE 
L = 1 

DO 16 J=1,NPHASE 
READ(10,*) QQ 
DO 17 K= L,L+NFAPH(J)-1 
A(I,K) = QQ 

17 CONTINUE 

L = L + NFAPH(J) 

16 CONTINUE 

15 CONTINUE 

ELSE 

C INPUT THE PROBABILITY OF SUCCESS IN A SINGLE TRIAL FOR EACH CAUSE 
C IN THE FIRST PHASE IF FRELOP EQUALS ZERO. 

DO 20 I=1,NCAUSE 

READ(10,*) A(I,1) 

20 CONTINUE 

ENDIF 

INPUT THE REMAINING VARIABLES , THE NUMBER OF SUCCESSFUL TRIALS 
BEFORE A DISCOUNT IS APPLIED (N); THE DISCOUNT FACTOR (R); THE SEED 
FOR THE RANDOM NUMBER GENERATOR, GGUBFS, (DSEED); RELIABILITY 
GROWTH FRACTION (FRIMP); TRIGGER FOR PRINTING INTERMEDIATE OUTPUT 
(lOPT) 

TRIGGERS FOR SAVING EACH ESTIMATE AT EACH PHASE FOR EACH ESTIMATOR 
lOPTl : MLE WITH DISCOUNTING 
I0PT2 : SINGLE PHASE MLE 
I0PT3 : EXPONENTIAL REGRESSION MODEL 

I0PT4 : WEIGHTED EXPONENTIAL REGRESSION MODEL (METHOD 1) 

I0PT5 : WEIGHTED EXPONENTIAL REGRESSION MODEL (METHOD 2) 

I0PT6 : WEIGHTED EXPONENTIAL REGRESSION MODEL (METHOD 3) 

DISCOUNTING OPTION TRIGGER (DISOPT); LLOYD FAILURE DISCOUNTING 

PARAMETER (GAMMA); LLOYD DISCOUNT INTERVAL 

READ(10,*) N 
READdO,*) R 
READ(10,*) DSEED 
READdO,*) FRIMP 
READdO,*) NREP 
READdO,*) lOPT 
READdO,*) lOPTl 
READdO,*) I0PT2 
READdO,*) I0PT3 
READdO,*) I0PT4 
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READCIO,*) IOPT5 
READCIO,*) IOPT6 
READCIO,*) DISOPT 
READCIO,*) GAMA 
READCIO,*) EDI 

C READ IN THE USER INPUTTED WEIGHTS FOR COMPUTATION OF METHOD 3. 

DO 25 I=1,NPHASE 

READCIO,*) WC3,I) 

25 CONTINUE 

XNREP = NREP 
DSEEDl = DSEED 

C INITIALIZE THE ARRAYS USED TO COMPUTE THE MEAN AND STANDARD DEVIATION 

C OF EACH ESTIMATOR 

DO 30 J=1,NPHASE 
DO 30 1=1,4 

MLEWDCI,J) 

MLESPCI,J) 

REGESTCI,J) 

WREST1CI,J) 

WREST2CI, J) 

WREST3CI,J) 

30 CONTINUE 



DO 31 J=1,NPHASE 
DO 31 1=1,6 

PHRESTCI,J) = 0. 0 

31 CONTINUE 



YJKCl) =0.0 
DO 35 1=1,999 

YJKCI+1) = YJKCl) + 1. O/I 
35 CONTINUE 

C COMPUTE AND STORE THE YSTAR VALUES UP TO 1000 
DO 39 J=1,NPHASE 

YSTARCJ,!) = l./NFAPHCJ) 

DO 41 I = 1,999 

YSTARCJ, I+l)= YSTARCJ,!) + 1. 0/CNFAPHC J)+I) 
41 CONTINUE 
39 CONTINUE 



C PRINT OUT THE YSTAR MATRIX FOR THE FIRST 219 ENTRIES OF THE lOTH PHASE 
WRITE C90,*) CYSTARC10,J), J=l,219) 



= 0.0 
= 0 . 0 
= 0 . 0 

= 0 . 0 

= 0 . 0 

= 0.0 



C COMPUTE AND STORE K BAR FOR THE EXPONENTIAL REGRESSION MODEL 



SUM = 0. 0 

DO 50 I=1,NPHASE 
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SUM = SUM + I 
REG(1,I) = SUM/I 
50 CONTINUE 

C MAJOR REPETITION OF THE SIMULATION LOOP 
DO 55 REP=1,NREP 
55 CONTINUE 



DO 500 REP=1,NREP 



C INITIALIZE FAILURE CAUSE VECTOR (NFCAUS) 

C COMPUTE THE INITIAL SYSTEM RELIABILITY 

REL = 1. 

DO 60 I=1,NCAUSE 
NFCAUS(I) = 0 
REL = REL * A(I,1) 

60 CONTINUE 

C INITIALIZE COLUMN (FAILURE # ) COUNTER FOR THE WORKING ARRAY (A) 
J = 1 

C LOOP TO COMPUTE THE NUMBER OF TRIALS UP TO AND INCLUDING FAILURE 
C AND THE CAUSE OF FAILURE FOR EACH FAILURE IN EACH PHASE 

DO 130 K=1,NPHASE 

C SKIP ACTUAL COMPONENT RELIABILITY COMPUTATION AFTER FIRST REP 
C AND FOR FIRST FAILURE 

IF(J. EQ. 1) GOTO 75 
IF(REP. GT. 1) GOTO 75 
REL = 1. 



C IF FIXED RELIABILITY OPTION IS SELECTED THEN PHASE RELIABILITIES 
C ARE COMPUTED AS FOLLOWS 

IF (FRELOP .EQ. 1) THEN 
DO 65 I=1,NCAUSE 
REL = REL*A(I,J) 

NFCAUS(I) = 0 
65 CONTINUE 

ELSE 

C COMPUTE NEW ACTUAL RELIABILITY FOR THE COMPONENT IN PHASE K 
DO 70 I=1,NCAUSE 

C INCREASE CAUSE PR( SUCCESS) IF IT CAUSED FAILURE IN THE PREVIOUS PHASE 
C COMPUTE NEXT PHASE RELIABILITY AND REINITIALIZE NFCAUS (NOT USED IF 
C FIXED PHASE RELIABILITY OPTION IS SELECTED). 



73 



IF(NFCAUS(I).EQ. 1) THEN 

A(I,J) = A(I,(J-D) + ((1. - A(I,(J-1)))*FRIMP) 
ELSEIF(NFCAUS(I). NE. 1) THEN 
A(I,J) = A(I,(J-D) 

ELSE 
END IF 

REL = REL*A(I,J) 

NFCAUS(I) = 0 
70 CONTINUE 



ENDIF 



75 J1 = 1 

TRTOT = 0. 0 

C COMPUTE THE NUMBER OF TRIALS UP TO AND INCLUDING FAILURE AND THE 

C CAUSE OF FAILURE FOR EACH FAILURE IN THE PHASE 

DO 120 L=1,NFAPH(K) 

IF(REP. GT. 1) GOTO 90 
IF(J1. EQ. 1) GOTO 85 
IF (FRELOP .EQ. 1) GOTO 85 
DO 80 I=1,NCAUSE 
A(I,J) = A(I,(J-D) 

80 CONTINUE 

85 A(ISYSPR,J) = REL 

A(IPHASE,J) = K 

90 MIN=7.2E75 

DO 110 I=1,NCAUSE 

C ASSIGN # TRIALS FOR CAUSES WITH PR( SUCCESS) = 0 OR 1 

IF(A(I,J).GE. 1. ) THEN 
A((I+NCAUSE),J) = 7.2E75 
GOTO 100 

ELSEIF(A(I,J). EQ. 0. ) THEN 
A((I+NCAUSE),J) = 1. 

GOTO 100 
ELSE 
ENDIF 



C CONVERT UNIFORM (0,1) RANDOM VARIABLE TO GEOMETRIC (// TRIALS UNTIL 
C FAILURE ) FOR EACH FAILURE CAUSE. RECORD THE MIN TRIALS FOR THE 
C CAUSES AS THE SYSTEM TRIALS UP TO AND INCLUDING FAILURE AND 
C RECORD THE FAILURE CAUSE 



A((I+NCAUSE),J) = INT(l.+(LOG(GGUBFS(DSEED))/LOG(A(I,J)))) 
100 IF(A((I+NCAUSE) ,J). LE. MIN) THEN 

MIN = A((I+NCAUSE),J) 

IMIN = I 
ELSE 
ENDIF 

110 CONTINUE 

A(IFAILC,J) = IMIN 
NFCAUS(IMIN) = 1 
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C COMPUTE THE TOTAL # OF TRIALS FOR THE MLE SINGLE PHASE ESTIMATE AND 
C INCREMENT FAILURE COUNTERS 

A(INTR,J) = MIN 

TRTOT = TRTOT + A(INTR,J) 

J = J + 1 
J1 = J1 + 1 
120 CONTINUE 

C COMPUTE THE MLE ESTIMATE OF COMPONENT RELIABILITY FOR THIS PHASE AND 
C COMPUTE THE RUNNING SUM OF ESTIMATES AND THE SUM OF ESTIMATES SQUARED 
C FOR COMPUTATION OF THE MEAN AND STANDARD DEVIATION OF THE ESTIMATE 

PHREST(2,K) = (TRTOT - NFAPH(K) ) /TRTOT 
MLESP(1,K) = MLESP(1,K) + PHREST(2,K) 

MLESP(2,K) = MLESP(2,K) + ( PHREST( 2 ,K)**2) 

130 CONTINUE 

C INITIALIZE THE ADJUSTED NUMBER OF FAILURES TO 1 AND THE COUNT OF THE 
C NUMBER OF TRIALS SINCE FAILURE OR DISCOUNTING (NTRIALS(I) ) TO 0 
C IN PREPARATION FOR THE DISCOUNTING ROUTINE 

DO 140 J=1,NFAIL 

A(IADJF,J) = 1. 

140 CONTINUE 

DO 150 I=1,NCAUSE 
NTRIAL(I) = 0 
150 CONTINUE 

C DISCOUNTING ROUTINE REVIEWS ALL PAST FAILURES AND CAUSES TO DATE 
C AND DETERMINES IF THE DISCOUNTING CONDITIONS HAVE BEEN MET. COMPUTES 
C THE ADJUSTED FAILURES, THE ADJUSTED // OF TRIALS AND YJK 

C INITIALIZE THE TOTAL ADJUSTED TRIALS IN A PHASE VECTOR, TADJTP; 

C INITIALIZE THE SUM OF THE VARIANCES,SVAR,THE VARIANCES , VAR 
C AND THE WEIGHTS, W FOR THE WEIGHTED EXPONENTIAL REGRESSION MODEL. 

DO 155 I = 1, NPHASE 

TADJTP (I) = 0 
SVARl(I) = 0 
SVAR2(I) = 0 
SUMW(I) = 0 
VAR(I) = 0 
W(1,I) = 0 
W(2,I) = 0 

155 CONTINUE 



SUMIA = 0 
J = 0 
LL = 1 

DO 300 K=l, NPHASE 
DO 200 L=1,NFAPH(K) 
J = J + 1 
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C UPDATES THE NUMBER OF TRIALS SINCE FAILURE OR DISCOUNTING FOR EACH 
C FAILURE CAUSE 

ICAUSE = INT(A(IFAILC,J)+. 5) 

DO 160 I=1,NCAUSE 

IF( ICAUSE. EQ. I) THEN 
NTRIAL(I) = 0 
ELSEIFC ICAUSE. NE. I) THEN 

NTRIAL(I) = NTRIAL(I) + INT(A( INTR, J)+. 5) 

ELSE 

ENDIF 

160 CONTINUE 



200 CONTINUE 

C CHOOSE DISCOUNTING METHOD TO BE USED 
IF(DISOPT. NE. 2) GOTO 180 
C PERFORM LLOYD'S FAILURE DISCOUNTING METHOD 



DO 170 1=1, J 

II = INT(A(IFAILC,I)+. 5) 

IF(NTRIAL(I1). EQ. 0) THEN 
A(IADJF,I) = 1. 0 
GOTO 170 

ELSE 

ENDIF 

C THIS IS THE MODIFIED LLOYD METHOD USING A DISCOUNT INTERVAL. THE 
C ORIGINAL DISCOUNT METHOD MAY BE EMPLOYED BY SETTING LDI TO O.NE. 

ALD = INT(NTRIAL(I1)/LDI) 

IF(ALD .EQ. 0) THEN 
A(IADJF,I) =1.0 
GO TO 170 
ELSE 

A(IADJF,I) = 1.0 - ((1. -GAMA)**(1. 0/ALD)) 

ENDIF 



170 CONTINUE 
GOTO 210 



C PERFORMS STRAIGHT PERCENT FAILURE DISCOUNTING AND 
C COMPUTES THE ADJUSTED OF FAILURES 

180 DO 190 1=1, J 

II = INT(A(IFAILC,I)+. 5) 

IF(NTRIAL(I1).EQ. 0) THEN 
A(IADJF,I) = 1. 

ELSEIF(NTRIAL(I1).GE.N) THEN 

A(IADJF,I) = A(IADJF,I)*((1. -R)**(NTRIAL(I1)/N)) 
ELSE 
ENDIF 

190 CONTINUE 
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C ADJUSTS THE # TRIALS SINCE FAILURE OR DISCOUNTING FOR THOSE CAUSES 
C THAT HAVE MET OR SURPASSED THE DISCOUNTING THRESHOLD 
C FOR THE STRAIGHT PERCENT DISCOUNTING METHOD 

DO 205 I=1,NCAUSE 

IF(NTRIALd).GE.N) NTRIAL(I) = MOD(NTRIAL( I) ,N) 

205 CONTINUE 

210 TADJT = 0. 0 

C COMPUTES THE ADJUSTED # OF TRIALS FROM THE ADJUSTED # OF FAILURES 
C AND COMPUTES THE SUM OF THE ADJUSTED # OF TRIALS FOR ESTIMATE COMP. 

C IF ADJUSTED FAILURES ARE APPROACHING 0 THEN ADJUSTED TRIALS MUST 
C BE PRE-SET. 

DO 240 I = 1,J 

IF(A(IADJF,I) .LE. .0000001) THEN 
A(IADJF,I) = .0000001 

END IF 

A(IADJT,I) = A(INTR,I)/A(IADJF,I) 

TADJT = TADJT + A(IADJT,I) 

240 CONTINUE 



C COMPUTE THE ADJUSTED NUMBER OF TRIALS IN A PHASE 

DO 245 M = LL,(LL+NFAPH(K) - 1) 

TADJTP(K) = TADJTP(K) + A(IADJT,M) 

245 CONTINUE 

LL = LL + NFAPH(K) 



C COMPUTE THE MLE ESTIMATE OF PHASE RELIABILITY USING DISCOUNTING 
PHREST(1,K) = (TADJT - J)/TADJT 

C COMPUTE Y-BAR AND Y-BAR FOR THE PHASE USING THE YSTAR MATRIX IF 
C THE ADJUSTED NUMBER OF TRIALS IS LESS THAN OR EQUAL TO 1000 
C AND USING EULER'S APPROXIMATION IF THE VALUE IS GREATER THAN 1000. 

IF ( NINT(TADJTP(K)) . EQ. NFAPH(K) )THEN 
REG(3,K) = 0.0 

ELSE IF (NINT(TADJTP(K) - NFAPH(K) ) . LE. 1000 ) THEN 
REG(3,K) = YSTAR(K, NINT(TADJTP(K) -NFAPH(K) ) ) 

ELSE 

IX = NINT(TADJTP(K)-1 ) 

X = IX 
Q = 12*X 
T = X +1 
S = X+2 

U = (EUL + (LOG(X)) + (1/(2*X)) - (1/(Q*T)) -(1/(Q*T*S)) ) 
REG(3,K) = U - YJK(NFAPH(K)-1) 

END IF 
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SUMIA = SUMIA + REG(3,K) 
REG(2,K) = SUMIA/K 



C COMPUTE THE EXPONENTIAL REGRESSION ESTIMATE BEGINNING WITH B HAT 



SUM = 0. 0 
SUMS = 0. 0 

IF (K. EQ. 1) GOTO 252 
DO 250 I = 1,K 

SUM = SUM + ((I-REG(1,K))*REG(3,I)) 

SUMS = SUMS + (d-REG(l,K))**2) 

250 CONTINUE 

REG(4,K) = SUM/SUMS 

C COMPUTE A HAT 

REG(5,K) = REG(2,K) - (REG( 4 ,K)*REG( 1 ,K) ) 

C COMPUTE AND STORE THE EXPONENTIAL REGRESSION ESTIMATE 

PHREST(3,K) = 1.0 - EXP( -(REG(5 ,K) + (REG( 4 ,K)*K) ) ) 
IF(PHREST(3,K). LT. 0. 0) PHREST(3,K)=0. 0 
GOTO 255 ' 

252 PHREST(3,K) = 1. 0 - EXP( -REG( 3 , 1) ) 

IF(PHREST(3,K). LT. 0. 0) PHREST( 3 ,K)=0. 0 



C STORE THE RUNNING SUM OF THE ESTIMATES FOR THE CURRENT PHASE AND THE 
C RUNNING SUN OF THE ESTIMATES SQUARED FOR COMPUTATION OF THE MEAN AND 
C STANDARD DEVIATION OF EACH ESTIMATE FOR EACH RELIABILITY GROWTH 
C MODEL 

255 MLEWD(1,K) = MLEWD(1,K) + PHREST(2,K) 

MLEWD(2,K) = MLEWD(2,K) + (PHREST(2,K)**2) 

REGEST(1,K) = REGEST(1,K) + PHREST(3,K) 

REGEST(2,K) = REGEST(2,K) + (PHREST(3,K)**2) 

C COMPUTE THE VARIANCES OF THE UNWEIGHTED EXPONENTIAL REGRESSION 
C ESTIMATES AND STORE THE RUNNING SUM OF THE VARIANCES FOR USE IN THE 
C WEIGHTED REGRESSION MODEL. 

VAR(K) = (NFAPH(K)*PHREST(1,K))/((1. -PHREST(1,K))**2) 

IF ( VAR(K) .LT. . 0000001)THEN 
VAR(K) = . 0000001 
END IF 

IF (K .EQ. 1) GO TO 258 

SVARl(K) = SVARl(K-l) + (l./VAR(K)) 

SVAR2(K) = SVAR2(K-1) + VAR(K) 

SWW(K) = SUMW(K-l) + W(3,K) 

GO TO 259 

258 SVARl(K) = (l./VAR(l)) 

SVAR2(K) = VAR( 1) 

SUM\v’(K) =W(3,1) 
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C COMPUTE THE WEIGHTS FOR EACH PHASE 
259 W(1,K) = (1. /VAR(K))/SVAR1(K) 

W(2,K) = VAR(K)/SVAR2(K) 

C COMPUTE AND STORE X-BAR(W) AND Y-BAR(W) 



SUMXl 
SUMYl 
SUMX2 
SUMY2 
SUMX3 
SUMY3 
DO 261 I 
SUMXl 
SUMYl 
SUMX2 
SUMY2 
SUMX3 
SUMY3 



1,K 

SUMXl 

SUMYl 

SUMX2 

SUMY 2 

SUMX3 

SUMY3 



+ (1. /VAR(K))*I 
+ (1. /VAR(K))*REG(3,I) 
+ VAR(K)*I 
+ VAR(K)*REG(3,I) 

+ W(3,K)*I 
+ W(3,K)*REG(3,I) 



261 CONTINUE 

WRECK 1,K) = SUMXl/SVARKK) 
WREG1(4,K) = SUMYl/SVARKK) 
WREG2(1,K) = SUMX2/SVAR2(K) 
WREG2(4,K) = SUMY2/SVAR2(K) 
WREG3(1,K) = SUMX3/SUMVs'(K) 
WREG3(4,K) = SUMY3/SUMW(K) 

C COMPUTE AND STORE B-HAT(W) 



SUMl =0.0 
SUMSl = 0.0 
SUM2 =0.0 
SUMS2 = 0. 0 
SUM3 =0.0 
SUMS3 = 0.0 

IF (K .EQ. 1) GO TO 272 
DO 271 I = 1,K 

SUMl = SUMl + (W(1,K)*(I-WREG1(1,K))*REG(3,I)) 
SUMSl = SUMSl + ( W(1,K)*((I-WREG1(1,K))**2)) 
SUM2 = SUM2 + (W(2,K)*(I-WREG2(1,K))*REG(3,I)) 
SUMS2 = SUMS2 + ( W( 2 ,K)*( ( I -WREG2( 1 ,K) )**2) ) 
SUM3 = SUM3 + (W(3,K)*(I-WREG3(1,K))*REG(3,I)) 
SUMS3 = SUMS3 + ( W( 3 ,K)*( ( I-WREG3( 1 ,K) )**2) ) 
271 CONTINUE 



WRECK 2, K) = SUMl/SUMSl 
WREG2(2,K) = SUM2/SUMS2 
WREG3(2,K) = SUM3/SUMS3 

C COMPUTE AND STORE A-HAT(W) 

WREG1(3,K) = WREGK4.K) - (WRECK 2 ,K)*WREG1( 1 ,K) ) 
WREG2(3,K) = WREG2(4,K) - ( WREG2( 2 ,K)*WREG2( 1 ,K) ) 
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WREG3(3,K) = WREG3(4,K) - (WREG3(2,K)*WREG3( 1,K)) 

C COMPUTE AND STORE THE WEIGHTED EXPONENTIAL REGRESSION ESTIMATE 

PHREST(4,K) = 1.0 - EXP( -(WREG1( 3 ,K)+(WREG1( 2 ,K)*K) ) ) 

PHREST(5,K) = 1.0 - EXP( -(WREG2( 3 ,K)+(WREG2( 2 ,K)*K) ) ) 

PHREST(6,K) = 1.0 - EXP( -(WREG3( 3 ,K)+(WREG3( 2,K)*K) ) ) 

IF (PHREST(4,K) . LT. 0.0) PHREST(4,K) =0.0 
IF (PHREST(5,K) . LT. 0.0) PHREST(5,K) = 0.0 
IF (PHREST(6,K) . LT. 0.0) PHREST(6,K) =0.0 
GO TO 275 

272 PHREST(4,K) = 1. 0 - EXP( -WREG1(4, 1) ) 

IF (PHREST(4,K) . LT. 0.0) PHREST(4,K) =0.0 

PHREST(5,K) = 1.0 - EXP( -WREG2(4, 1) ) 

IF (PHREST(5,K) . LT. 0.0) PHREST(5,K) =0.0 

PHREST(6,K) = 1.0 - EXP( -WREG3(4, 1) ) 

IF (PHREST(6,K) . LT. 0.0) PHREST(6,K) = 0.0 

C STORE THE RUNNING SUM OF THE WEIGHTED ESTIMATES FOR THE CURRENT 
C PHASE AND THE RUNNING SUM OF THE ESTIMATES SQUARED FOR COMPUTATION OF 
C THE MEAN AND STANDARD DEVIATION OF EACH WEIGHTED ESTIMATE FOR THE 
C WEIGHTED EXPONENTIAL REGRESSION MODEL. 
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WRESTl(l.K) = WRESTl(l.K) + PHREST(4,K) 
WREST1(2,K) = WREST1(2,K) + (PHREST(4,K) 



2 ) 



WREST2(1,K) = WREST2(1,K) + PHREST(5,K) 
WREST2(2,K) = WREST2(2,K) + ( PHREST(5 ,K)**2) 



WREST3(1,K) = WREST3(1,K) + PHREST(6,K) 
WREST3(2,K) = WREST3(2,K) + (PHREST( 6 ,K)**2) 



C STORE THE ACTUAL PHASE RELIABILITY 



AREL(K) = A(ISYSPR,J) 

C PRINT INTERMEDIATE OUTPUT IF REQUESTED AND THE NUMBER OF REPETITIONS 
C IS NOT GREATER THAN 5 

IFCIOPT. NE. 1) GOTO 300 
IF(REP.GT. 5) GOTO 300 

WRITE(25,1000' REP,K 

1000 FORMAT(T16, 'REPETITION NUMBER; ',14,' PHASE NUMBER: ',14) 
WRITE(25,1010) A(ISYSPR,J) 

1010 FORMATC 2 2X, 'ACTUAL COMPONENT RELIABILITY: ',F7.5) 

WRITE(25,1022) PHREST(1,K) 

1022 FORMATC 2 OX, 'MLE ESTIMATE USING DISCOUNTING; ',F7.5) 

WRITE(25,1025) PHREST(2,K) 

1025 FORMATC 18X,' MLE ESTIMATE OF PHASE RELIABILITY: ',F7.5) 

WRITEC25,1027) PHRESTC4,K) 

1027 FORMATC 14X, 'WEIGHTED REG. ESTIMATE C METHOD 1 ) ; ',F7.5) 



SO 



WRITE(25,1026) PHREST(5,K) 

1026 F0RMAT(14X, 'WEIGHTED REG. ESTIMATE ( METHOD 2 ) : ’,F7.5) 

WRITE(25,1029) PHREST(6,K) 

1029 F0RMAT(14X, 'WEIGHTED REG. ESTIMATE ( METHOD 3 ) : ',F7.5) 

WRITE(25,1028) PHREST(3,K) 

1028 F0RMAT(14X, 'REGRESSION ESTIMATE OF PHASE RELIABILITY: ',F7.5) 
WRITE(25,1030) 

1030 FORMATC ' ' , ' ' ) 

DO 260 I=1,NCAUSE 

WRITE(25 , 1035)1 ,A( I , J) ,A( ( I+NCAUSE) , J) 

1035 FORMATC 12X, 'CAUSE: ,13,' PR(SUCCESS): ',F7.6,' //TRIALS: ', 

CFIO. 0) 

260 CONTINUE 

WRITE(25,1036) 

1036 FORMATC ' ' , ' ' ) 

WRITEC25,1040) 

1040 FORMATC 4X, 'FAIL //',3X,'FAIL CAUSE ', 3X ,' // TRIALS' ,3X, 'AD J // FAIL' ,3 
CX,'ADJ // TRIALS') 

DO 270 1=1, J 

WRITEC25,1050)I,ACIFAILC,I) ,ACINTR,I) ,ACIADJF,I) ,ACIADJT,I) 

1050 F0RMATC4X,I3,8X,F3. 0 , 7X,F8. 0 ,4X,F8. 6 ,4X, F12. 0,3X,F11. 4) 

270 CONTINUE 

WRITEC25,1060) 

1060 FORMATC ' ',///) 

300 CONTINUE 



C PRINT EACH OF THE ESTIMATES TO THEIR APPROPRIATE OUTPUT FILE 
C IF REQUESTED 



401 

402 

403 

404 

405 

406 
2000 



IFCIOPTl.NE. 1) 
WRITEC55,2000) 
IFCI0PT2.NE. 1) 
WRITEC65.2000) 
IFCI0PT3.NE. 1) 
WRITEC75,2000) 
IFCI0PT4. NE. 1) 
WRITEC 15,2000) 
IFCI0PT5.NE. 1) 
WRITEC39,2000) 
IFCI0PT4. NE. 1) 
WRITEC49,2000) 
FORMATC ' ' , 30C 



GOTO 402 
CPHRESTCl,!) , 
GOTO 403 
CPHRESTC2,!) , 
GOTO 404 
CPHRESTC3,!) , 
GOTO 405 
CPHRESTC4,!) , 
GOTO 406 
CPHRESTC5,!) , 
GOTO 500 
CPHRESTC6,!) , 
'7. 6: IX)) 



I=1,NPHASE) 

I=1,NPHASE) 

I=1,NPHASE) 

I=1,NPHASE) 

I=1,NPHASE) 

I=1,NPHASE) 



500 CONTINUE 



C PRINT OUT THE WORKING "A" MATRIX IN MATRIX LISTING A1 
DO 4050 J = l,C2*NCAUSE)+6 
WRITEC81,*) CACJ,I),I = 1,NFAIL) 

4050 CONTINUE 



C PRINT OUT THE THE TOTAL NUMBER OF TRIALS TO SYSTEM FAILURE IN EACH 
C PHASE 

DO 4051 J = l.NPHASE 

WRITEC89,*) 'T ',TADJTPCJ) 
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4051 CONTINUE 



C PRINT OUT THE WEIGHTS COMPUTED BY METHOD 1 TO THE TRIALS FILE 
DO 4052 J = 1,NPHASE 

WRITE(89,*) 'Ml ' ,W(1,J) 

4052 CONTINUE 



C PRINT OUT THE WEIGHTS COMPUTED BY METHOD 2 TO THE TRIALS FILE 
DO 4053 J = 1,NPHASE 

WRITE(89,*) 'M2 ' ,W(2,J) 

4053 CONTINUE 



C PRINT OUT THE REGRESSION PARAMETERS IN MATRIX FORM TO THE REGMAT FILE 



DO 4054 J = 1,5 

WRITE(87,*) (REG(J,I),I = 1,NPHASE) 
4054 CONTINUE 



C PRINT OUT THE WEIGHTED REGRESSION PARAMETERS IN MATRIX FORM TO THE 
C REGMAT FILE 

DO 4055 J = 1,4 

WRITE(87,*) (WREG1(J,I),I = 1,NPHASE) 

4055 CONTINUE 

DO 4056 J = 1,4 

WRITE(87,*) (WREG2(J,I),I = 1,NPHASE) 

4056 CONTINUE 



DO 4057 J = 1,4 

WRITE(87,*) (WREG3(J,I),I = 1,NPHASE) 

4057 CONTINUE 

C UPON COMPLETION OF ALL REPETITIONS, COMPUTE THE MEAN AND STANDARD 
C DEVIATION OF EACH ESTIMATE FOR EACH PHASE SKIPPING COMPUTATIONS IF 
C ONLY ONE REPETITION IS REQUIRED 

IF (NREP. LE. 1) GOTO 601 

DO 600 I=1,NPHASE 
MLEWD(3,I) = MLEWD(1,I)/XNREP 
MLESP(3,I) = MLESP(1,I)/XNREP 
REGEST(3,I) = REGEST(1,I)/XNREP 
WREST1(3,I) = WREST1(1,I)/XNREP 
WREST2(3,I) = WREST2(1,I)/XNREP 
WREST3(3,I) = WREST3(1,I)/XNREP 

MLEWD( 4 , I ) = SQRT( ( MLEWD( 2 , I ) - ( XNREP*( MLEWD( 3 , I ) **2 ) ) ) / ( XNREP - 1 ) ) 
MLESP(4,I) = SQRT((MLESP(2,I)-(XNREP*(MLESP(3,I)**2)))/(XNREP-1)) 
REGE ST( 4 , 1 ) =SQRT( ( REGEST( 2 , 1 ) - ( XNREP*( REGE ST( 3 , 1 ) **2 ) ) ) / ( XNRE P - 1 ) ) 
WRESTK4, 1)=SQRT( (WRESTK 2,1) -(XNREP*(WREST1( 3 , 1)**2) ) )/(XNREP- 1) ) 
WREST2(4,I)=SQRT((WREST2(2,I)-(XNREP*(WREST2(3,I)**2)))/(XNREP-1)) 
WREST3 ( 4 , I ) =SQRT( ( WRE ST3 ( 2 , I ) - ( XNREP*( WRE ST3 ( 3 , I) **2 ) ) ) / ( XNREP - 1 ) ) 
600 CONTINUE 



C PRINT THE FINAL OUTPUT TABLE TO A FILE 



601 WRITE(35,3000) 
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3000 FORMATC 'o' ,T47, 'DISCRETE RELIABILITY GROWTH SIMULATION') 
WRITE(35,3010) 

3010 FORMATC '-' ,T5 4, 'MODEL PARAMETER SUMMARY') 

WRITE(35,3020) NCAUSE 

3020 FORMATC 'o' ,T4 7, 'NUMBER OF POSSIBLE FAILURE CAUSES ',14) 

IF CFRELOP .EQ. 1) GOTO 4000 
WRITEC35,3030) 

3030 FORMATC 'O' ,T38, 'CAUSE NUMBER' ,T64, ' SINGLE TRIAL PRC SUCCESS ) FOR 
CPHASE 1') 

DO 3050 M=l, NCAUSE 
WRITEC35,3040) M,ACM,1) 

3040 FORMATC ' ' ,T43 , 12 ,T79 ,F8. 6) 

3050 CONTINUE 

WRITEC35,3060) FRIMP 

3060 FORMATC 'O' ,T37, 'FRACTION CAUSE RELIABILITY IMPROVES AFTER FAILURE 
C' ,F8. 6) 

5000 WRITEC35,3080) NPHASE 

3080 FORMATC'-' ,T48,' NUMBER OF PHASES IN THE SIMULATION ',12) 
WRITEC35,3090) 

3090 FORMATC 'o' ,T42, 'PHASE NUMBER' ,T59 , 'NUMBER OF FAILURES IN THE FIRST 
C PHASE') 

DO 3110 M=l, NPHASE 
WRITEC35,3100) M,NFAPHCM) 

3100 FORMATC' ' ,T43 , 12 ,T73 , 12) 

3110 CONTINUE 

WRITEC35,3120) NFAIL 

3120 FORMATC ' O' ,T51,' TOTAL NUMBER OF FAILURES ',14) 

IFCDISOPT. EQ. 2) GO TO 3160 
WRITEC35,3130) 

3130 FORMATC '-' ,T38, 'DISCOUNTING PERFORMED USING THE CONSTANT FRACTION 
CMETHOD' ) 

WRITEC35,3140) R 

3140 FORMATC 'O' ,T44, 'FRACTION EACH FAILURE IS DISCOUNTED ' ,F8. 6) 
WRITEC35,3150) N 

3150 FORMATC' ' ,T33 ,' NUMBER OF TRIALS AFTER A FAILURE BEFORE A DISCOUNT 
C IS APPLIED ' ,14) 

GO TO 3190 
3160 WRITEC35.3170) 

3170 FORMATC'-' ,T44, 'DISCOUNTING PERFORMED USING THE LLOYD METHOD') 
WRITEC35,3180) GAMA 

3180 FORMATC 'O' ,T39,' PERCENT C. I. ( USED AS DISCOUNT FRACTION ) ' ,F8. 6 
C) 

WRITEC35,3185) LDI 

3185 FORMATC 'o' ,T50, 'LLOYD DISCOUNT INTERVAL: ' ,13) 

3190 WRITEC35,3200) DSEEDl 

3200 FORMATC '-' ,T46, 'RANDOM NUMBER SEED USED ' ,F15. 2) 

WRITEC35,3210) NREP 

3210 FORMATC 'o' ,T3 7, 'NUMBER OF REPETITIONS OF THE SIMULATION PERFORMED 
C',I7) 

WRITEC35,3220) 

3220 FORMATC ' l' ,T61, 'ESTIMATOR: ' ) 

WRITEC35,3230) 

3230 FORMATC 'o' ,T48, 'SINGLE PHASE MLE WITHOUT DISCOUNTING') 
WRITEC35,3240) 

3240 FORMATC'-' ,T60, 'MEAN' ,T83,' ESTIMATE ' ,T109,' 95 %' ) 

WRITEC35,3250) 
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3250 FORMATC ’,T12, 'phase NUMBER' ,T29 ACTUAL RELIABILITY' ,T52 ,' PREDIC 
CTED RELIABILITY' ,T78, 'STANDARD DEVIATION' ,T101 ,' CONFIDENCE INTERVA 
CL' ) 

C COMPUTE C. I. FOR SINGLE PHASE MLE 



DO 3270 M=1,NFHASE 

Cl = (1. 96*MLESP(4,M))/SQRT(XNREP) 

CIU = MLESP(3,M) + Cl 
CIL = MLESP(3,M) - Cl 

WRITE(35,3260)’ M,AREL(M) ,MLESP( 3 ,M) ,MLESP( 4 ,M) , CIL,CIU 
3260 FORMATC 'o' ,T17,I2,T34,F8. 6,T58,F8. 6,T82,F9. 6,T99, '( ',F8.6,' , ' ,F 
C8.6,' )') 

3270 CONTINUE 

WRITE(35,3220) 

WRITE(35,3280) 

3280 FORMATC 'o' ,T42, 'MAX LIKELIHOOD ESTIMATE USING DISCOUNTED FAILURES' 
C) 

WRITEC35,3240) 

WRITEC35,3250) 

C COMPUTE C. I. FOR MLE WITH DISCOUNTING 



DO 3290 M=1,NPHASE 

Cl = Cl- 96*MLEWDC4,M))/SQRTCXNREP) 

CIU = MLEWDC3.M) + Cl 
CIL = MLEWDC3.M) - Cl 

WRITEC35.3260) M,ARELCM) ,MLEWDC3,M) ,MLEWDC4,M) ,CIL,CIU 
3290 CONTINUE 



WRITEC35.3220) 

WRITEC35,3320) 

3320 FORMATC 'O' ,T43, 'REGRESSION ESTIMATE USING DISCOUNTED FAILURES') 
WRITEC35,3240) 

WRITEC35,3250) 

C COMPUTE C. I. FOR EXPONENTIAL REGRESSION ESTIMATES 



DO 3330 M=1,NPHASE 

Cl = C 1- 96*REGESTC4,M))/SQRTCXNREP) 

CIU = REGESTC3,M) + Cl 
CIL = REGESTC3,M) - Cl 

WRITEC35,3260) M.ARELCM) ,REGESTC3,M) ,REGESTC4,M) ,CIL,CIU 
3330 CONTINUE 

C WEIGHTED REGRESSION CMETHOD 1) 

WRITEC35,3220) 

WRITEC35,3321) 

3321 FORMATC 'O' ,T4 3, 'WEIGHTED REGRESSSION ESTIMATE CMETHOD 1 )') 

WRITEC35,3240) 

WRITEC35,3250) 

C COMPUTE C. I. FOR WEIGHTED EXPONENTIAL REGRESSION ESTIMATES 



DO 3331 M=1,NPHASE 
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Cl = (1. 96*WREST1(4,M))/SQRT(XNREP) 

CIU = WREST1(3,M) + Cl 
CIL = WREST1(3,M) - Cl 

WRITE(35,3260) M,AREL(M) ,WREST1(3,M) ,WREST1(4,M) ,CIL,CIU 

3331 CONTINUE 

C WEIGHTED REGRESSION (METHOD 2) 

WRITE(35,3220) 

WRITE(35,3322) 

3322 FORMATC ' O’ ,T4 3, 'WEIGHTED REGRESSION ESTIMATE (METHOD 2) ’) 

WRITE(35,3240) 

WRITE(35,3250) 

C COMPUTE C. I. FOR WEIGHTED EXPONENTIAL REGRESSION ESTIMATES 

DO 3332 M=1,NPHASE 

Cl = (1. 96*WRi:ST2(4,M))/SQRT(XNREP) 

CIU = WREST2(3,M) + Cl 
CIL = WREST2(3,M) - Cl 

WRITE(35,3260) M,AREL(M) ,WREST2(3,M) ,WREST2(4,M) , CIL, CIU 

3332 CONTINUE 

C WEIGHTED REGRESSION (METHOD 3) 

WRITE(35,3220) 

WRITE(35,3323) 

3323 FORMAT( 'o’ ,T4 3, 'WEIGHTED REGRESSION ESTIMATE (METHOD 3) ') 

WRITE(35,3240) 

WRITE(35,3250) 

C COMPUTE C.I. FOR WEIGHTED EXPONENTIAL REGRESSION ESTIMATES 



DO 3333 M=1,NPHASE 

Cl = (1. 96*WREST3(4,M))/SQRT(XNREP) 

CIU = WREST3(3,M) + Cl 
CIL = WREST3(3,M) - Cl 

WRITE(35,3260) M,AREL(M) ,WREST3(3,M) ,WREST3(4,M) ,CIL,CIU 
3333 CONTINUE 



WRITE(35,3340) 

3340 FORMAT( ' l’ ,T59 , 'RECAPITULATION' //) 

WRITE(35,3350) 

3350 FORMAT( ,T3, 'PHASE' Til, 'ACTUAL' ,T28,'MEAN' ,T38,'EST' ,T53, 
C’MEAN' ,T63,'EST' ,T78,‘MEAN' ,T88,'EST' ,T103,'MEAN’ ,T113,‘EST^) 
WRITE(35,3360) 

3360 FORMAT( ' ’ ,T11 'RELIAB' ,T28, 'WT' ,T38 ,'STD' ,T53, 'MLE' ,T63, 'STD' ,T7 
C7, 'PHASE' ,T88,'STD' ,T103,'REG' ,T113,'STD') 

WRITE(35,3370) 

3370 FORMATC ' ,T28,'REG' ,T35 ,' DEVIATION’ ,T53,'W/D' ,T60 ,’ DEVIATION' 
C,T78, 'MLE' ,T85, 'DEVIATION' ,T103, 'EST' ,T110, 'DEVIATION' ) 
WRITE(35,3375) 

3375 FORMATC ' ,T28 , ' EST’ / ) 

DO 650 I=1,NPHASE 

WRITE(35,3380) I ,AREL( I ) ,WREST1( 3 , I ) ,WREST1(4 , I ) ,MLEWD( 3 , I ) 

C , MLEWD( 4,1), MLESP( 3,1), MLESP( 4,1) ,REGEST( 3,1), REGEST( 4,1) 
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3380 F0RMAT( 'O' ,T4 , 13 ,T11 ,F7. 6 ,T26 ,F7. 6 ,T36 ,F7. 6,T51,F7. 6,T61,F7. 6 
C,T76,F7. 6,T86,F7. 6,T101,F7. 6,T111,F7. 6) 

650 CONTINUE 

C PAGE 2 OF RECAPITULATION: 

WRITE(35,3381) 

3381 FORMATC ' 1' ,T59, 'RECAPITULATION CONT. '//) 

VRITE(35,3382) 

3382 FORMATC ' ,T3, 'PHASE* .Til , 'ACTUAL* ,T28, 'MEAN* ,T38, 'EST' .T53, 
C'MEAN' ,T63,'EST' ,T78,'MEAN' ,T88,'EST' ,T103,'MEAN' ,T113,'EST') 

WRITE(35,3383) 

3383 FORMATC* ' .Til, 'RELIAB' ,T28 'REG* ,T38, ' STD* ,T53, 'METHOD* ,T63, 
C'STD' ,T77 , 'method' ,T88, ' STD ' ,T103, 'METHOD* ,T113,'STD') 

WRITEC35,3384) 

3384 FORMATC ’ ' ,T28,'EST' ,T35 ,' DEVIATION* ,T53,'0NE* ,T60 ,' DEVIATION* 
C,T78,’TW0* ,T85, 'DEVIATION* ,T103 THREE ' ,T110 DEVIATION* ) 

DO 651 I=1,NPHASE 

WRITEC 35 , 3385 ) I , ARELC I ) , REGESTC 3,1) ,REGESTC 4,1), WRESTlC 3,1) 

C , WRESTIC 4 I ) , WREST2C 3,1) , WREST2C 4 , 1 ) ,WREST3C 3,1) ,WREST3C 4,1) 

3385 FORMATC *0' ,T4 , 13 ,T1 1 ,F7. 6,T26 ,F7. 6 ,T36 ,F7. 6 ,T5 1 , F7. 6,T6l,F7. 6 
C,T76,F7. 6,T86,F7. 6,T101,F7. 6,T111,F7. 6) 



C PRINT RELIABILITY ESTIMATES TO FILES FOR PLOTTING PURPOSES 

WRITEC51,3400) MLEWDC4,I) 

WRITEC 5 2, 3400) MLEWDC3,I) 

WRITEC 61, 3400) MLESPC4,!) 

WRITEC 62, 3400) MLESPC3,!) 

WRITEC 7 1,3400) REGESTC 4, I) 

WRITEC 72, 3400) REGESTC 3, I) 

WRITEC 76,3400) WRESTIC 4, I) 

WRITEC 77,3400) WREST1C3,I) 

WRITEC 78,3400) WREST2C4,!) 

WRITEC 79, 3400) WREST2C3,I) 

WRITEC91,3400) WREST3C4,!) 

WRITEC 92, 3400) WREST3C3,I) 

WRITEC 82 3400) ARELC I) 

3400 FORMATC ' ,F7. 6) 

651 CONTINUE 



GO TO 6000 
4000 WRITEC 35, 4010) 

4010 FORMATC IX, //,T50,' FIXED PHASE RELIABILITY OPTION*) 
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WRITE(35,4020) 

4020 FORMAT( ,T38, 'PHASE NUMBER' ,T7 8, 'ACTUAL RELIABILITY') 
DO 4030 H=1,NPHASE 

WRITE(35,4040) M,AREL(M) 

4040 FORMATC 'O' ,T41,I2,T83,F8. 6) 

4030 CONTINUE 
GO TO 5000 
6000 CONTINUE 

STOP 

END 
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4. JIMC FORTRAN PROGRAM 



Vc Vc V% Vc V? V? Vr Vr Vc Vv Vf V? Vc Vc Vc Vr Vc Vc Vf Vc Vv Vf Vc Vc Vc /V Vc V? V* V% Vf Vc Vr Vc 



* 



DISCRETE RELIABILITY GROVsTH SIMULATION * 

* 

PROGRAMMED BY JAMES E DRAKE, * 

JAMES D CHANDLER, * 

AND PAM A MARKIEWICZ 

LAST MODIFIED 11 JUN 1988 * 



THE FOLLOWING EXTERNAL FILES ARE USED BY THE PROGRAM * 

INPUT : DATA AND PARAMETER INPUT FILE (DEVICE # 10) * 

THESIS : OUTPUT FILE CONTAINING INTERMEDIATE COMPUTATIONS * 

(DEVICE /) 20) * 

RELIAB: OUTPUT FILE CONTAINING FINAL RESULTS OF THE SIMULATION * 

(DEVICE /) 30) * 

EST : OUTPUT FILE CONTAINING EACH PHASE ESTIMATE FOR EACH * 

REPLICATION OF THE WOODS WEIGHTED AVERAGE ESTIMATE * 

(DEVICE ()40) * 



MLEWD : OUTPUT FILE CONTAINING MLE ESTIMATES USING DISCOUNTING * 
FOR EACH PHASE AND EACH REPLICATION * 

(DEVICE 50) * 

MLESP : OUTPUT FILE CONTAINING MLE ESTIMATE FOR EACH SINGLE PHASE * 
AND ALL REPLICATIONS USING NO DISCOUNTING 

(DEVICE 60) * 

REGEST ; OUTPUT FILE CONTAINING EACH PHASE ESTIMATE FOR EACH * 

REPLICATION OF THE EXPONENTIAL REGRESSION ESTIMATE * 

(DEVICE # 70) * 

THE FOLLOWING IS A LIST OF KEY ARRAYS USED IN THE SIMULATION 



A : MAIN WORKING ARRAY CONTAINS PROBABILITY OF SUCCESS FOR 

EACH FAILURE CAUSE, NUMBER OF TRIALS UNTIL FAILURE FOR * 
EACH FAILURE CAUSE AND THE SYSTEM, CAUSE OF FAILURE, 

PHASE NUMBER, ADJUSTED NUMBER OF TRIALS AND ADJUSTED 
NUMBER OF FAILURES * 

DIMENSION ( ((2*//CAUSES)+6),()FAILURES ) 

NFAPH : CONTAINS THE NUMBER OF FAILURES IN EACH PHASE * 

DIMENSION (1,/) PHASES) * 

NFCAUS : BINARY ARRAY USED TO DETERMINE IF A FAILURE OCCURRED IN * 
A PHASE * 

DIMENSION ( 1, (/FAILURE CAUSES) * 

NTRIAL : CONTAINS THE NUMBER OF TRIALS SINCE LAST FAILURE OR * 

DISCOUNTING FOR EACH FAILURE CAUSE * 

DIMENSION ( 1, (/FAILURE CAUSES ) * 

TADJTP : CONTAINS THE TOTAL NUMBER OF ADJUSTED TRIALS IN A PHASE * 
DIMENSION ( 1, (/ OF PHASES) * 

PHREST : RECORDS THE PHASE ESTIMATE FOR EACH ESTIMATOR WITHIN A * 
SINGLE REPLICATION * 

DIMENSION (7, (/PHASES) * 

ROW 1 : WOODS WEIGHTED AVERAGE EST * 

ROW 2 : MLE WITH DISCOUNTING * 
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* 




ROW 3 


SINGLE PHASE MLE 


* 


Vf 




ROW 4 


WOODS REGRESSION ESTIMATE (UNWEIGHTED) 


* 


* 




ROW 5 


WEIGHTED REGRESSION ESTIMATE (METHOD 1) 


* 


* 




ROW 6 


WEIGHTED REGRESSION ESTIMATE (METHOD 2) 


* 


* 




ROW 7 


WEIGHTED REGRESSION ESTIMATE (METHOD 3) 


•it 


Vr 


AREL : 


CONTAINS ACTUAL SYSTEM RELIABILITY IN EACH PHASE 


* 


* 




DIMENSION (1, ^/PHASES) 


* 


* 


CUMSF : 


CONTAINS THE NUMBER OF SUCCESS AND FAILURES FOR EACH 


•it 


* 




FAILURE CAUSE (USED WITH WOODS WEIGHTED AVERAGE EST. ) 


•it 


* 




DIMENSION (B./^FAILUE CAUSES) 


* 






ROW 1 


NUMBER OF FAILURES 


* 


* 




ROW 2 


NUMBER OF SUCCESSES 


* 


* 




ROW 3 


ADJUSTED NUMBER OF SUCCESSES 


* 


Vc 


YJK : 


CONTAINS YJK VALUES UP TO 1000 


* 


* 




DIMENSION (1,1000) 


Vc 


* 


VAR : 


CONTAINS THE VARIANCE OF THE ESTIMATE AT EACH PHASE 


Vc 






DIMENSION (1, //PHASES) 


* 


* 


SVARl : 


CONTAINS THE RUNNING SUM OF 1/VAR(K) UP TO THE KTH PHASE 




* 




DIMENSION (1, //PHASES) 


* 


* 


SVAR2 : 


CONTAINS THE RUNNING SUM OF VAR(K) UP TO THE KTH PHASE 


* 


* 




DIMENSION (1, //PHASES) 


* 


* 


SUMW : 


CONTAINS THE RUNNING SUM OF WEIGHTS UP TO THE KTH PHASE 


* 


* 




DIMENSION (1, //PHASES) 


* 


Vf 


W : 


CONTAINS THE WEIGHT GIVEN TO EACH RELIABILITY ESTIMATE AT 


Vc 


Vc 




EACH PHASE. DIMENSION (3, //PHASES) 


Vc 


* 




ROW 1: 


WEIGHT = (1/VAR(K))/(SUM OF 1/VAR(K) UP TO 


Vc 


* 






THE KTH PHASE) 


Vc 


* 




ROW 2; 


WEIGHT = VAR(K)/(SUM OF VAR(K) UP TO KTH PHASE) 


Vc 


* 




ROW 3: 


WEIGHT = USER INPUTTED VALUES. 


Vc 


* 


REG : 


ARRAY USED TO COMPUTE THE EXPONENTIAL REGRESSION ESTIMATE 


Vc 


* 




DIMENSION (5, //PHASES) 


* 






ROW 1 


K BAR 


Vc 


/V 




ROW 2 


Y BAR 


Vc 


* 




ROW 3 


Y BAR FOR THE PHASE 


Vc 


* 




ROW 4 


B HAT 


Vc 


Vc- 




ROW 5 


A HAT 




Vf 


WREGl : 


ARRAY USED TO COMPUTE THE WEIGHTED EXPONENTIAL REGRESSION 


•it 






ESTIMATE USING METHOD 1 WEIGHTS. DIMENSION (4, //PHASES) 


Vc 


Vf 




ROW 1 


X BAR(W) 


•it 


J* 




ROW 2 


B HAT(W) 


V" 


•it 




ROW 3 


A HAT(W) 


Vc 


•it 




ROW 4 


Y BAR(W) 


* 


•it 


WREG2 : 


SAME AS ABOVE USING METHOD 2 WEIGHTS 


Vc 


•it 


WREG3 ; 


SAME AS ABOVE USING METHOD 3 WEIGHTS 


•it 




THE REMAINING ARRAYS ARE USED TO COMPUTE THE MEAN AND VARIANCE 


•it 


•it 


OF EACH 


ESTIMATE AT EACH PHASE. THEY ALL HAVE THE SAME DIMENSIONS 


•it 


•it 


AND STRUCTURE 




V" 


•it 




DIMENSION (4, //PHASES) 


* 


•it 




ROW 1 


RUNNING SUM OF ESTIMATES 


Vc 


•it 




ROW 2 


RUNNING SUM OF SQUARED ESTIMATES 


•it 


•it 




ROW 3 


MEAN OF THE ESTIMATES 


•it 


•it 

"it 




ROW 4 


STANDARD DEVIATION OF THE ESTIMATES 


Vc 


•it 


MLEWD : 


VALUES 


FOR THE MLE WITH DISCOUNTING 


•it 


mX. 


MLESP : 


VALUES 


FOR THE SINGLE PHASE MLE 


•it 
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* REGEST : VALUES FOR THE EXPONENTIAL REGRESSION ESTIMATE * 

* WRESTl : VALUES FOR THE WEIGHTED EXPONENTIAL REGRESSION ESTIMATE 

* USING METHOD 1 WEIGHTS. * 

* WREST2 : VALUES FOR THE WEIGHTED EXPONENTIAL REGRESSION ESTIMATE 

* USING METHOD 2 WEIGHTS. * 

* WRESTS : VALUES FOR THE WEIGHTED EXPONENTIAL REGRESSION ESTIMATE * 

* USING METHOD 3 WEIGHTS. * 






C DEFINE AND DIMENSION VARIABLES 



PARAMETER (NR=50 ,NC=200) 

INTEGER REP,DISOPT,FRELOP,LDI,ALD 
REAL*4 MIN 

REALMS DSEED , MLE SP , MLEWD , EUL , SUMl , SUMS 1 , SUM2 , SUMS2 , SUMS 
C,SUMS3 

DIMENSION NFAPH ( NR ) , A( NR , NO) , NFC AUS ( NR ) , NTRI AL ( NR ) , PHRE ST ( 7 , NR ) , 
CMLEWD(4,NR) ,MLESP(4,NR) ,REGEST(4,NR) ,AREL(NR) ,YJK( 1000) , 
CREG(5,NR) ,TADJTP(1000) , VAR(NR) , SVARl(NR) ,W(3,NR) , 

CWREG1(4,NR) ,WREST1(4,NR) ,WREST2(4,NR) ,WREST3(4,NR) , 
CSVAR2(NR),WREG2(4,NR) ,WREG3(4,NR) ,SUMW(NR) ,EST(4,NR) ,CUMSF(3,NR) 

C READ IN THE NUMBER OF CAUSES TO BE USED ( NCAUSE ) AND THE NUMBER 
C OF PHASES ( NPHASE ) IN THE TEST 

RE AD (10,*) NCAUSE 
RE AD (10,*) NPHASE 

C CHECK IF FIXED RELIABILITY OPTION IS DESIRED. FIX EULER'S NUMBER. 



READ(10,*) FRELOP 
EUL = 0. 5772156648 

C CREATE VARIABLES FOR THE ROW INDICES OF THE WORKING MATRIX ( A ) 

C I PHASE: PHASE 

C ISYSPR: ACTUAL C0MP0NT;NT RELIABILITY 

C INTR: NUMBER OF TRIALS UP TO AND INCLUDING FAILURE 
C IFAILC: CAUSE OF THE FAILURE 

C lADJF: ADJUSTED NUMBER OF FAILURES ED 

C AFTER DISCOUNTING HAS BEEN APPLIED 

C lADJT: ADJUSTED NUMBER OF TRIALS AFTER DISCOUNTING HAS BEEN APPLIED 

C lYJK: YJK COMPUTED ON THE ADJUSTED NUMBER OF TRIALS 

I PHASE = (2*NCAUSE)+1 
ISYSPR = IPHASE +1 
INTR = ISYSPR + 1 
IFAILC = INTR + 1 
lADJF = IFAILC + 1 
lADJT = lADJF + 1 
lYJK = lADJT + 1 
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nnonnonnnoononn 



C READ IN THE NUMBER OF FAILURES IN EACH PHASE ( NFAPH(I) ) AND 
C COMPUTE THE TOTAL NUMBER OF FAILURES IN THE TEST ( NFAIL ) 

NFAIL = 0 

DO 10 I=1,NPHASE 

READ(10,*) NFAPH(I) 

NFAIL = NFAIL + NFAPH(I) 

10 CONTINUE 

C INPUT THE PROBABILITY OF SUCCESS IN A SINGLE TRIAL FOR EACH CAUSE 
C IN EACH PHASE IF FRELOP EQUALS ONE. 

IF (FRELOP .EQ. 1) THEN 
DO 15 I=1,NCAUSE 
L =1 

DO 16 J=1,NPHASE 
READ(10,*) QQ 

DO 17 K = L,L+NFAPH(J)-1 
A(I,K) = QQ 

17 CONTINUE 

L = L+NFAPH(J) 

16 CONTINUE 

15 CONTINUE 

ELSE 

C INPUT THE PROBABILITY OF SUCCESS IN A SINGLE TRIAL FOR EACH CAUSE 
C IN THE FIRST PHASE IF FRELOP EQUALS ZERO. 

DO 20 I=1,NCAUSE 

READCIO,*) A(I,1) 

20 CONTINUE 



END IF 

INPUT THE REMAINING VARIABLES , THE NUMBER OF SUCCESSFUL TRIALS 
BEFORE A DISCOUNT IS APPLIED (N); THE DISCOUNT FACTOR (R); THE SEED 
FOR THE RANDOM NUMBER GENERATOR, GGUBFS , (DSEED); RELIABILITY 
GROWTH FRACTION (FRIMP); TRIGGER FOR PRINTING INTERMEDIATE OUTPUT 
(lOPT) 

TRIGGERS FOR SAVING EACH ESTIMATE AT EACH PHASE FOR EACH ESTIMATOR 
lOPTl : WOODS WEIGHTED AVERAGE MODEL 
I0PT2 : MLE WITH DISCOUNTING 
I OPTS : SINGLE PHASE MLE 
I0PT4 : EXPONENTIAL REGRESSION MODEL 

I0PT5 : WEIGHTED EXPONENTIAL REGRESSION MODEL (METHOD 1) 

I0PT6 ; WEIGHTED EXPONENTIAL REGRESSION MODEL (METHOD 2) 

I0PT7 : WEIGHTED EXPONENTIAL REGRESSION MODEL (METHOD 3) 

DISCOUNTING OPTION TRIGGER (DISOPT); LLOYD FAILURE DISCOUNTING 
PARAMETER (GAMMA); LLOYD DISCOUNT INTERVAL 

READdO,*) N 
READdO,*) R 
READdO,*) DSEED 
READdO,*) FRIMP 
READdO,*) NREP 
READdO,*) lOPT 
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READ(10,*) lOPTl 
READ(10,*) IOPT2 
READCIO,*) IOPT3 
READ(10,*) IOPT4 
READ(10,*) IOPT5 
READ(10,*) IOPT6 
READ(10,*) IOPT7 
READ(10,*) DISOPT 
READ (10,*) GAMA 
READ(10,*) LDJ 

C READ IN THE USER INPUTTED WEIGHTS FOR COMPUTATION OF METHOD 3. 

DO 25 I=1,NPHASE 

READ(10,*) W(3,I) 

25 CONTINUE 

XNREP = NREP 
DSEEDl = DSEED 

C INITIALIZE THE ARRAYS USED TO COMPUTE THE MEAN AND STANDARD DEVIATION 

C OF EACH ESTIMATOR 

DO 30 J=1,NPHASE 
DO 30 1=1,4 

EST(I,J) = 0. 0 
MLEWD(I,J) = 0. 0 
MLESP(I,J) = 0. 0 
REGEST(I,J) = 0. 0 
WREST1(I,J) = 0.0 
WREST2(I,J) =0.0 
WREST3(I,J) =0.0 

30 CONTINUE 

DO 31 J=1,NPHASE 
DO 31 1=1,6 

PHREST(I,J) = 0. 0 

31 CONTINUE 



C COMPUTE AND STORE THE YJK VALUES UP TO 1000 

YJK(l) = 0. 0 
DO 40 1=1,999 

YJK(I+1) = YJK(I) + 1. O/I 
40 CONTINUE 

C COMPUTE AND STORE K BAR FOR THE EXPONENTIAL REGRESSION MODEL 

SUM = 0. 0 
DO 50 I=1,NPHASE 
SUM = SUM + I 
REG(1,I) = SUM/ I 
50 CONTINUE 

C MAJOR REPETITION OF THE SIMULATION LOOP 
DO 500 REP=1,NREP 
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C INITIALIZE FAILURE CAUSE VECTOR (NFCAUS) AND (CUMSF) 
C COMPUTE THE INITIAL SYSTEM RELIABILITY 



REL = 1. 

DO 60 I=1,NCAUSE 
NFCAUSCI) = 0 
REL = REL * A(I,1) 

DO 60 J=l,3 

CUMSF(J,I) = 0 

60 CONTIKTJE 

C INITIALIZE COLUMN (FAILURE # ) COUNTER FOR THE WORKING ARRAY (A) 

J = 1 

C LOOP TO COMPUTE THE NUMBER OF TRIALS UP TO AND INCLUDING FAILURE 
C AND THE CAUSE OF FAILURE FOR EACH FAILURE IN EACH PHASE 

DO 130 K=1,NPHASE 

C SKIP ACTUAL COMPONENT RELIABILITY COMPUTATION AFTER FIRST REP 
C AND FOR FIRST FAILURE 

IF(J. EQ. 1) GOTO 75 
IF(REP. GT. 1) GOTO 75 
REL = 1. 

C IF FIXED RELIABILITY OPTION IS SELECTED THEN PHASE RELIABILITIES 
C ARE COMPUTED AS FOLLOWS 

IF (FRELOP .EQ. 1) THEN 
DO 65 I=1,NCAUSE 
REL = REL*A(I,J) 

NFCAUS(I) = 0 
65 CONTINUE 

ELSE 

C COMPUTE NEW ACTUAL RELIABILITY FOR THE COMPONENT IN PHASE K 
DO 70 I=1,NCAUSE 

C INCREASE CAUSE PR( SUCCESS) IF IT CAUSED FAILURE IN THE PREVIOUS PHASE 
C COMPUTE NEXT PHASE RELIABILITY AND REINITIALIZE NFCAUS (NOT USED IF 
C FIXED PHASE RELIABILITY OPTION IS SELECTED). 

IF(NFCAUS(I).EQ. 1) THEN 

A(I,J) = A(I,(J-1)) + ((1. - A(I,(J-1)))*FRIMP) 
ELSEIF(NFCAUSd). NE. 1) THEN 
A(I,J) = A(I,(J-1)) 

ELSE 

ENDIF 

REL = REL*A(I,J) 

NFCAUS(I) = 0 
70 CONTINUE 

ENDIF 
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75 J1 = 1 

TRTOT = 0.0 

C COMPUTE THE NUMBER OF TRIALS UP TO AND INCLUDING FAILURE ANT) THE 
C CAUSE OF FAILURE FOR EACH FAILURE IN THE PHASE 

DO 120 L=1,NFAPH(K) 

IF(REP. GT. 1) GOTO 90 
IF(J1. EQ. 1) GOTO 85 
IF (FRELOP .EQ. 1) GOTO 85 
DO 80 I=1,NCAUSE 
A(I,J) = A(I,(J-1)) 

80 CONTINUE 

85 A(ISYSPR,J) = REL 

A(IPHASE,J) = K 

90 MIN = 7. 2E75 

DO 110 I=1,NCAUSE 

C ASSIGN TRIALS FOR CAUSES WITH PR( SUCCESS) = 0 OR 1 

IF(A(I,J).GE. 1. ) THEN 
A((I+NCAUSE),J) = 7.2E75 
GOTO 100 

ELSEIF(A(I,J). EQ. 0. ) THEN 
A((I+NCAUSE),J) = 1. 

GOTO 100 
ELSE 
END IF 

C CONVERT UNIFORM (0,1) RANDOM VARIABLE TO GEOMETRIC (// TRIALS UNTIL 
C FAILURE ) FOR EACH FAILURE CAUSE. RECORD THE MIN TRIALS FOR THE 
C CAUSES AS THE SYSTEM fr TRIALS UP TO AND INCLUDING FAILURE AND 
C RECORD THE FAILURE CAUSE 

A((I+NCAUSE),J) = INT(l.+(LOG(GGUBFS(DSEED))/LOG(A(I,J)))) 
100 IF(A((I+NCAUSE),J).LE.MIN) THEN 

MIN = A((I+NCAUSE) , J) 

IMIN = I 
ELSE 
END IF 

110 CONTINUE 

A(IFAILC,J) = IMIN 
NFCAUS(IMIN) = 1 

C COMPUTE THE TOTAL # OF TRIALS FOR THE MLE SINGLE PHASE ESTIMATE AND 
C INCREMENT FAILURE # COUNTERS 

A(INTR,J) = MIN 

TRTOT = TRTOT + A(INTR,J) 

J = J + 1 
J1 = J1 + 1 
120 CONTINUE 



94 



C COMPUTE THE MLE ESTIMATE OF COMPONENT RELIABILITY FOR THIS PHASE AND 
C COMPUTE THE RUNNING SUN OF ESTIMATES AND THE SUM OF ESTIMATES SQUARED 
C FOR COMPUTATION OF THE MEAN AND STANDARD DEVIATION OF THE ESTIMATE 

PHREST(3,K) = (TRTOT - NFAPH(K) ) /TRTOT 
MLESP(1,K) = MLESP(1,K) + PHREST(3,K) 

MLESP(2,K) = MLESP(2,K) + ( PHREST( 3 , K)**2) 

130 CONTINUE 

C INITIALIZE THE ADJUSTED NUMBER OF FAILURES TO 1 AND THE COUNT OF THE 
C NUMBER OF TRIALS SINCE FAILURE OR DISCOUNTING (NTRIALS(I) ) TO 0 
C IN PREPARATION FOR THE DISCOUNTING ROUTINE 

DO 140 J=1,NFAIL 

A(IADJF,J) = 1. 

140 CONTINUE 

DO 150 I=1,NCAUSE 
NTRIAL(I) = 0 
150 CONTINUE 

C INITIALIZE THE VARIANCE, VAR, THE SUM OF 1 / VAR, S VAR, AND THE WEIGHTS, W. 

DO 155 I=1,NPHASE 

SVARl(I) = 0 
SVAR2(I) = 0 
SUMW(I) = 0 
VAR(I) = 0 
W(1,I) = 0 
W(2,I) = 0 

155 CONTINUE 

C DISCOUNTING ROUTINE REVIEWS ALL PAST FAILURES AND CAUSES TO DATE 
C AND DETERMINES IF THE DISCOUNTING CONDITIONS HAVE BEEN MET. COMPUTES 
C THE ADJUSTED FAILURES, THE ADJUSTED OF TRIALS AND YJK 

J = 0 

DO 300 K=1,NPHASE 
DO 200 L=1,NFAPH(K) 

J = J + 1 

C UPDATES THE NUMBER OF TRIALS SINCE FAILURE OR DISCOUNTING FOR EACH 
C FAILURE CAUSE 

ICAUSE = INT(A(IFAILC,J)+. 5) 

DO 160 I=1,NCAUSE 

IF( ICAUSE. EQ. I) THEN 
NTRIAL(I) = 0 
ELSEIFCICAUSE. NE. I) THEN 

NTRIAL(I) = NTRIAL(I) + INT( A( INTR, J)+. 5) 

ELSE 
END IF 

160 CONTINUE 
200 CONTINUE 
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C CHOOSE DISCOUNTING METHOD TO BE USED 



IF(DISOPT.NE. 2) GOTO 180 
C PERFORM LLOYD'S FAILURE DISCOUNTING METHOD 



DO 170 1=1, J 

II = INT(A(IFAILC,I)+. 5) 

IF(NTRIAL(I1). EQ. 0) THEN 
A(IADJF,I) = 1.0 
GOTO 170 

ELSE 
END IF 

C THIS IS THE MODIFIED LLOYD METHOD USING A DISCOUNT INTERVAL. THE 
C ORIGINAL DISCOUNT METHOD MAY BE EMPLOYED BY SETTING LDI TO ONE. 

ALD = INT(NTRIAL(I1)/LDI) 

IF(ALD .EQ. 0) THEN 
A(IADJF,I) = 1.0 
GO TO 170 
ELSE 

A(IADJF,I) = 1. 0 - ((1. -GAMA)**(1- 0/ALD)) 

ENDIF 



170 CONTINUE 
GOTO 210 



C PERFORMS STRAIGHT PERCENT FAILURE DISCOUNTING AND 
C COMPUTES THE ADJUSTED # OF FAILURES 

180 DO 190 1=1, J 

II = INT(A(IFAILC,I)+. 5) 

IF(NTRIAL(I1).EQ. 0) THEN 
A(IADJF,I) = 1. 

ELSEIF(NTRIAL(I1). GE.N) THEN 

A(IADJF,I) = A(IADJF,I)*(( 1. -R)**(NTRIAL(I1)/N)) 

ELSE 

ENDIF 

190 CONTINUE 

C ADJUSTS THE # TRIALS SINCE FAILURE OR DISCOUNTING FOR THOSE CAUSES 
C THAT HAVE MET OR SURPASSED THE DISCOUNTING THRESHOLD 
C FOR THE STRAIGHT PERCENT DISCOUNTING METHOD 

DO 205 I=1,NCAUSE 

IF(NTRIALd). GE.N) NTRIAL(I) = MOD(NTRIAL( I) ,N) 

205 CONTINUE 

210 TADJT = 0. 0 
TYJK = 0. 0 
TPYJK = 0. 0 
K1 = 0 



DO 215 12=1,3 
DO 215 I=1,NCAUSE 
CUMSF(I2,I) = 0 
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215 



CONTINUE 



C COMPUTES THE ADJUSTED # OF TRIALS FROM THE ADJUSTED OF FAILURES 
C AND COMPUTES THE SUM OF THE ADJUSTED # OF TRIALS FOR ESTIMATE COMP. 

PREL = 0. 0 
LTRIAL = 0 

C IF ADJUSTED FAILURES ARE APPROACHING 0 THEN ADJUSTED TRIALS MUST 
C BE PRE-SET. 

DO 240 1=1, J 

IF(A(IADJF,I) .LE. .0000001) THEN 
A(IADJF,I) = .0000001 
END IF 

A(IADJT,I) = A(INTR,I)/A(IADJF,I) 

TADJT = TADJT + A(IADJT,I) 

C COMPUTE YJK FROM THE ADJUSTED # OF TRIALS AND STORE THE SUM FOR 
C ESTIMATE COMPUTATION, USE ARRAY FOR TRIALS < 1000 AND APPROX. FOR 
C VALUES > 1000 

N1 = NINT(A(IADJT,I)) 

IF(N1. LE. 1000) THEN 

A(IYJK,I) = YJK(Nl) 

ELSEIF(N1. GT. 1000) THEN 
X=N1 
Q=12*X 
T=X+1 
S=X+2 

A( lYJK, I)=(EUL+( LOG(X) )+( l/( 2*X) ) -( 1/(Q*T) ) -( l/( Q*T*S) ) ) 

ELSE 

ENDIF 

C DETERMINE IF A PHASE BOUNDARY HAS BEEN REACHED TO BEGIN ESTIMATE 
C COMPUTATION 

IF(I. EQ. 1) GOTO 225 

IF( A( I PHASE , I ) . NE. A( IPHASE , ( I - 1 ) ) ) THEN 

C COMPUTE THE WOODS WEIGHTED AVERAGE ESTIMATE 

MAX = 0 
K1 = K1 + 1 

C DETERMINE THE FAILURE CAUSE WITH THE LARGEST # OF FAILURES 

DO 220 I1=1,NCAUSE 

IF(CUMSF( 1,11). GT. MAX) THEN 
MAX = CUMSF(1,I1) 

ICOL = II 
ELSE 
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END IF 
220 CONTINUE 

C COMPUTE YJK VALUE FOR THE CURRENT PHASE ESTIMATE 

IF(CUMSF(1,IC0L). LE. 1000) THEN 
AHATL = yJK(CUMSF(l,ICOL)) 

ELSEIF(CUMSF(l,ICOL). GT. 1000) THEN 
X = CUMSF(l,ICOL) 

Q=12*X 

T=X+1 

S=X+2 

AHATL=(EUL+(LOG(X))+(l/(2*X))-(l/(Q*T))-(l/(Q*'F'S))) 

ELSE 

ENDIF 

IX = CUMSF( l,ICOL) + CUMSF(3,ICOL) 

IF(IX. LE. 1000) THEN 
AHATU = YJK(IX) 

ELSEIFdX. GT. 1000) THEN 
X = IX 
Q=12*X 
T=X+1 
S=X+2 



AHATU=( EUL+( LOG(X) )+( l/( 2*X) ) - ( 1/ ( Q*T) ) - ( 1/ ( Q*'P''S ) ) ) 

ELSE 

ENDIF 



C COMPUTE CURRENT PHASE RELIABILITY ESTIMATE 

AHAT = AHATU - AHATL 
CREL = 1.0 - EXP(-AHAT) 

X = CUMSF(1,IC0L) + CUMSF(3,IC0L) 

C COMPUTE AND STORE THE WOODS WEIGHTED AVERAGE ESTIMATE 

PREL = ((LTRIAL*PREL)/X) + ( ( (X-LTRIAL)*CREL)/X) 

LTRIAL = CUMSF(1,IC0L) + CUMSF( 3 , ICOL) 

C COMPUTE THE PHASE AND GLOBAL AVERAGE FOR YJK USED IN THE EXPONENTIAL 
C REGRESSION ESTIMATES ARE 

REG(2,K1) = TYJK/d-l) 

REG(3,K1) = TPYJK/NFAPH(K1) 

TPYJK = 0. 0 



ENDIF 



C COMPUTE THE NUMBER OF FAILURES AND SUCCESSES FOR EACH FAILURE CAUSE 
C USED IN THE WOODS WEIGHTED AVERAGE ESTIMATE 

225 ICAUSE = INT(A(IFAILC,I)+. 5) 

DO 230 I1=1,NCAUSE 
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CUMSF(2,I1) = CUMSF(2,I1) + INT(A(INTR,I) + .5) 

CUMSF(3,I1'! = CUMSF(3,I1) + N1 
230 CONTINUE 

CUMSF(1,ICAUSE) = CUMSF( 1 , ICAUSE) + 1 
CUMSF(2,ICAUSE') = CUMSF( 2 , ICAUSE) - 1 
CUMSF( 3, ICAUSE) = CUMSF( 3 , ICAUSE) - 1 
TPYJK = TPYJK + A(IYJK,I) 

TYJK = TYJK + A(IYJK,I) 

240 CONTINUE 

C REPEAT COMPUTATIONS FOR THE WOODS WEIGHTED AVERAGE ESTIMATE FOR THE 
C FINAL PHASE 

MAX = 0 

K1 = K1 + 1 

DO 245 I1=1,NCAUSE 

IF(CUMSF(1,I1). GT. MAX) THEN 
MAX = CUMSF(1,I1) 

ICOL =11 
ELSE 
ENDIF 
245 CONTINUE 

IF(CUMSF( 1, ICOL). LE. 1000) THEN 
AHATL = YJK(CUMSF(l,ICOL)) 

ELSEIF(CUMSF(l,ICOL). GT. 1000) THEN 
X = CUMSF(l,ICOL) 

0=1 2*X 

T=X+1 

S=X+2 

AHATL=( EUL+C LOG( X ) ) +( 1 / ( 2*X ) ) - ( 1 / ( Q*T) ) - ( 1 / ( Q*T*S ) ) ) 

ELSE 

ENDIF 

IX = CUMSF(l,ICOL) + CUMSF(3,ICOL) 

IF(IX. LE. 1000) THEN 

AHATU = YJK(IX) 

ELSEIFCIX. GT. 1000) THEN 
X = IX 
0=1 2*X 
T=X+1 
S=X+2 

AHATU=( EUL+( LOG( X ) ) +( 1/ ( 2*X ) ) - ( 1 / ( 0*T) ) - ( 1 / ( 0*T’'^S ) ) ) 

ELSE 

ENDIF 

AHAT = AHATU - AHATL 
CREL = 1. 0 - EXP(-AHAT) 

X = CUMSFd.ICOL) + CUMSF(3,ICOL) 



PREL = ((LTRIAL*PREL)/X) + ( ( (X-LTRIAL)*CREL)/X) 
LTRIAL = CUMSF(l,ICOL) + CUMSF( 3 , ICOL) 
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REG(2,K1) = TYJK/(J) 

REG(3,K1) = TPYJK/NFAPH(K1) 

PHREST(1,K) = PREL 

C COMPUTE THE MLE ESTIMATE OF PHASE RELIABILITY USING DISCOUNTING 
PHREST(2,K) = (TADJT - J)/TADJT 

C COMPUTE THE EXPONENTIAL REGRESSION ESTIMATE BEGINNING WITH B HAT 

SUM = 0. 0 
SUMS = 0.0 

IF (K.EQ. 1) GOTO 252 
DO 250 I = 1,K 

SUM = SUM + ((I-REG(1,K))*REG(3,D) 

SUMS = SUMS + (d-REG(l,K))**2) 

250 CONTINUE 

REG(4,K) = SUH/SUMS 

C COMPUTE A HAT 

REG(5,K) = REG(2,K) - ( REG(4,K)*REG( 1 ,K) ) 

C COMPUTE AND STORE THE EXPONENTIAL REGRESSION ESTIMATE 

PHREST(4,K) = 1.0 - EXP( -(REG(5 ,K) + (REG(4,K)*K))) 
IF(PHREST(4,K). LT. 0. 0) PHREST(4,K)=0. 0 
GOTO 255 

252 PHREST(4,K) = 1. 0 - EXP( -REG(3, 1)) 

IF(PHREST(4,K). LT. 0. 0) PHREST(4,K)=0. 0 



C STORE THE RUNNING SUM OF THE ESTIMATES FOR THE CURRENT PHASE AND THE 
C RUNNING SUM OF THE ESTIMATES SQUARED FOR COMPUTATION OF THE MEAN AND 
C STANDARD DEVIATION OF EACH ESTIMATE FOR EACH RELIABILITY GROWTH 
C MODEL 

255 EST(1,K) = EST(1,K) + PHREST(1,K) 

EST(2,K) = EST(2,K) + (PHREST( 1 ,K)**2) 

MLEWD(1,K) = MLEWD(1,K) + PHREST(2,K) 

MLEWD(2,K) = MLEWD(2,K) + ( PHREST( 2 ,K)**2) 

REGEST(1,K) = REGEST(1,K) + PHREST(4,K) 

REGEST(2,K) = REGEST(2,K) + (PHREST(4,K)**2) 

C COMPUTE THE VARIANCES OF THE UNW’EIGHTED EXPONENTIAL REGRESSION 
C ESTIMATES AND STORE THE RUNNING SUM OF THE VARIANCES FOR USE IN THE 
C WEIGHTED REGRESSION MODEL. 

VAR(K) = (NFAPH(K)*PHREST(1,K))/((1. -PHREST(1,K))**2) 

IF ( VAR(K) .LT. . 0000001)THEN 
VAR(K) = .0000001 
END IF 

IF (K .EQ. 1) GO TO 258 

SVARl(K) = SVARl(K-l) + (l./VAR(K)) 
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SVAR2(K) = SVAR2(K-1) + VAR(K) 
SUMW(K) = SUMW(K-l) + W(3,K) 

GO TO 259 

258 SVARl(K) = (1. /VAR(l)) 

SVAR2(K) = VAR(l) 

SUMW(K) =W(3,1) 



C COMPUTE THE WEIGHTS FOR EACH PHASE 
259 W(1,K) = (1. /VAR(K))/SVAR1(K) 

W(2,K) = VAR(K)/SVAR2(K) 

C COMPUTE AND STORE X-BAR(W) AND Y-BAR(W) 
SUMXl = 0 
SUMYl = 0 
SUMX2 = 0 



SUMY2 = 0 
SUMX3 = 0 
SUMY3 = 0 
DO 261 I = 
SUMXl = 
SUMYl = 
SUMX2 = 
SUMY2 = 
SUMX3 = 
SUMY3 = 



l.K 

SUMXl + (1. /VAR(K))*I 
SUMYl + (1. /VAR(K))*REG(3,I) 
SUMX2 + VAR(K)*I 
SUMY2 + VAR(K)*REG(3,I) 

SUMX3 + V(3,K)*I 
SUHY3 + W(3,K)*REG(3,I) 



261 CONTINUE 

WRECK 1,K) = SUMXl/SVARKK) 
V;REG1(4,K) = SUMYl/SVARKK) 
WREG2(1,K) = SUMX2/SVAR2(K) 
WREG2(4,K) = SUMY2/SVAR2(K) 
WREG3(1,K) = SUMX3/SUMW(K) 
WREG3(4,K) = SUMY3/SUMWCK) 

C COMPUTE AND STORE B-HAT(W) 

SUMl =0.0 
SUMSl = 0.0 
SUM2 =0.0 
SUMS2 = 0.0 
SUM3 = 0. 0 
SUMS3 =0.0 

IF (K .EQ. 1) GO TO 272 



DO 271 I = 1,K 

SUMl = SUMl + (W(1,K)*(I-WREG1(1,K))*REG(3,I)) 
SUMSl = SUMSl + ( W(1,K)*((I-WREG1(1,K))**2)) 
SUM2 = SUM2 + (W(2,K)*(I-WREG2(1,K))*REG(3,I)) 
SUMS2 = SUMS2 + ( W( 2 ,K)*( ( I -WREG2( 1 ,K) )**2) ) 
SUM3 = SUM3 + (W(3,K)*(I-WREG3(1,K))*REG(3,D) 
SUMS3 = SUMS3 + ( W( 3 ,K)*( ( I -WREG3( 1 ,K) )**2) ) 
271 CONTINUE 



WREG1(2,K) = SUMl/SUMSl 
WREG2(2,K) = SUM2/SUMS2 
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WREG3(2,K) = SUM3/SUMS3 

C COMPUTE AND STORE A-HAT(W) 

WREG1(3,K) = WREG1(4,K) - (WREG1( 2 ,K)*WREG1( 1 ,K) ) 
WREG2(3,K) = WREG2(4,K) - (WREG2( 2 ,K)*WREG2( 1 ,K) ) 
WREG3(3,K) = WREG3(4,K) - (WREG3( 2 ,K)*WREG3( 1 ,K) ) 



C COMPUTE AND STORE THE WEIGHTED EXPONENTIAL REGRESSION ESTIMATE 

PHREST(5,K) = 1.0 - EXP( -(WREG1( 3,K)+(WREG1( 2 ,K)*K) ) ) 

PHREST(6,K) = 1.0 - EXP( -(WREG2( 3 ,K)+(WREG2( 2 ,K)*K) ) ) 

PHREST(7,K) = 1.0 - EXP( -(WREG3( 3,K)+( WREG3( 2 ,K)*K) ) ) 

IF (PHREST(5,K) . LT. 0.0) PHREST(5,K) =0.0 

IF (PHREST(6,K) . LT. 0.0) PHREST(6,K) =0.0 

IF (PHREST(7,K) . LT. 0.0) PHREST(7,K) = 0.0 

GO TO 275 

272 PHREST(5,K) = 1. 0 - EXP( -WREG1(4, 1) ) 

IF (PHREST(5,K) . LT. 0.0) PHREST(5,K) =0.0 

PHREST(6,K) = 1. 0 - EXP( -WREG2( 4 , 1) ) 

IF (PHREST(6,K) . LT. 0.0) PHREST(6,K) =0.0 

PHREST(7,K) = 1.0 - EXP( -WREG3( 4, 1) ) 

IF (PHREST(7,K) . LT. 0.0) PHREST(7,K) =0.0 

C STORE THE RUNNING SUM OF THE WEIGHTED ESTIMATES FOR THE CURRENT 
C PHASE AND THE RUNNING SUM OF THE ESTIMATES SQUARED FOR COMPUTATION OF 
C THE MEAN AND STANDARD DEVIATION OF EACH WEIGHTED ESTIMATE FOR THE 
C WEIGHTED EXPONENTIAL REGRESSION MODEL. 



275 WREST1(1,K) = WREST1(1,K) + PHREST(5,K) 

WREST1(2,K) = WREST1(2,K) + ( PHREST( 5 ,K)**2) 

WREST2(1,K) = WREST2(1,K) + PHREST(6,K) 
WREST2(2,K) = WREST2(2,K) + ( PHREST( 6 ,K)**2) 

WREST3(1,K) = WREST3(1,K) + PHREST(7,K) 
WREST3(2,K) = WREST3(2,K) + (PHREST( 7 ,K)**2) 



C STORE THE ACTUAL PHASE RELIABILITY 
AREL(K) = A(ISYSPR,J) 

C PRINT INTERMEDIATE OUTPUT IF REQUESTED AND THE NUMBER OF REPETITIONS 
C IS NOT GREATER THAN 5 

IFCIOPT. NE. 1) GOTO 300 
IF(REP. GT. 5) GOTO 300 

WRITE(20,1000) REP,K 
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1000 F0RMAT(T16,’ REPETITION NUMBER: ’,14,’ PHASE NUMBER: ’,14) 

WRITE(20,1010) A(ISYSPR,J) 

1010 FORMATC 2 2X, ’ACTUAL COMPONENT RELIABILITY: ’,F7.5) 

WRITE(20,1020) PHREST(1,K) 

1020 FORMATC 20X, 'PREDICTED COMPONENT RELIABILITY: ’,F7.5) 

WRITE(20,1022) PHREST(2,K) 

1022 FORMATC 20X, ’MLE ESTIMATE USING DISCOUNTING: ’,F7.5) 

WRITEC20,1025) PHRESTC3,K) 

1025 FORMATC 18X, ’MLE ESTIMATE OF PHASE RELIABILITY: ’,F7.5) 
WRITEC20,1027) PHRESTC4,K) 

1027 FORMATC 14X, ’REGRESSION ESTIMATE OF PHASE RELIABILITY: ’,F7.5) 

WRITEC20,1028) PHRESTC4,K) 

1028 FORMATC 14X, ’WEIGHTED REG. ESTIMATE C METHOD 1 ) : ’,F7.5) 

WRITEC20,1026) PHRESTC5,K) 

1026 FORMATC 14X, ’WEIGHTED REG. ESTIMATE C METHOD 2 ) : ’,F7.5) 

WRITEC20,1029) PHRESTC6,K) 

1029 FORMATC 14X, ’WEIGHTED REG. ESTIMATE C METHOD 3 ) : ’,F7.5) 

WRITEC20,1030) 

1030 FORMATC ’ ’ , ’ ' ) 

DO 260 I=1,NCAUSE 

WRITEC20,1035)I,ACI,J),ACCI+NCAUSE) ,J) 

1035 FORMATC 12X, ’CAUSE: ’,13,' PRCSUCCESS): ’,F7.6,’ //TRIALS: ’, 
CFIO. 0) 

260 CONTINUE 

WRITEC20,1036) 

1036 FORMATC ’ ’ , ’ ’ ) 

WRITEC20,1040) 

1040 FORMATC 4X, 'FAIL //’,3X,’FAIL CAUSE ’, 3X ,’ // TRIALS ’, 3X, ’ AD J // FAIL’, 3 
CX,'ADJ # TRIALS' ,7X, ’YJK’ ) 

DO 270 1=1, J 

WRITEC20,1050)I,ACIFAILC,I),ACINTR,I),ACIADJF,I) ,ACIADJT,I),ACIYJK 
C,I) 

1050 FORMATC 4X, 13, 8X, F3. 0,7X,F8. 0,4X,F8. 6,4X,F12. 0,3X,F11. 4) 

270 CONTINUE 

WRITEC20,1060) 

1060 FORMATC’ ’,///) 



300 CONTINUE 

C PRINT EACH OF THE 3 ESTIMATES TO THEIR APPROPRIATE OUTPUT FILE 
C IF REQUESTED 



400 

401 

402 

403 

404 

405 

406 
2000 



IFCIOPTl. NE. 1) GOTO 401 
WRITEC40,2000) CPHRESTCl,!) 
IFCIOPT2. NE. 1) GOTO 402 
WRITEC50,2000) CPHRESTC2,!) 
IFCIOPT3. NE. 1) GOTO 403 
WRITEC60,2000) C PHRESTC 3 , I ) , 
IFCI0PT4. NE. 1) GOTO 404 
WRITEC70,2000) CPHRESTC4,I) 
IFCIOPT5.NT;. 1) GOTO 405 
WRITEC16,2000) C PHRESTC 5, I) 
IFCIOPT5. NE. 1) GOTO 406 
WRITEC39,2000) C PHRESTC 6, I) 
IFCIOPT5.NE. 1) GOTO 500 
WRITEC49,2000) CPHRESTC7,!) 
FORMATC’ ’ ,30CF7. 6: IX)) 



I=1,NPHASE) 

I=1,NPHASE) 

I=1,NPHASE) 

I=1,NPHASE) 

I=1,NPHASE) 

I=1,NPHASE) 

I=1,NPHASE) 
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500 CONTINUE 



C WORKING A MATRIX 

DO 4050 J = l,(2*NCAUSE)+7 

WRITE(83,*) (A(J,I), I = 1,NFAIL) 
4050 CONTINUE 



DO 4052 J = 1,NCAUSE 

WRITE(88,*) (REG(J,I), I = 1,NPHASE) 

4052 CONTINUE 

C UPON COMPLETION OF ALL REPETITIONS, COMPUTE THE MEAN AND STANDARD 
C DEVIATION OF EACH ESTIMATE FOR EACH PHASE SKIPPING COMPUTATIONS IF 
C ONLY ONE REPETITION IS REQUIRED 

IF (NREP. LE. 1) GOTO 601 

DO 600 I=1,NPHASE 
EST(3,I) = EST(1,I)/XNREP 
MLEWD(3,I) = MLEVD(1,I)/XNREP 
MLESP(3,I) = MLESP(1,I)/XNREP 
REGEST(3,I) = REGEST(1,I)/XNREP 
WREST1(3,I) = WREST1(1,I)/XNREP 
WREST2(3,I) = WREST2C 1,I)/XNREP 
WREST3(3,I) = WREST3(1,I)/XNREP 

EST(4,I) = SQRT((EST(2,I)-(XNREP*(EST(3,I)**2)))/(XNREP-1)) 

MLEWD( 4 , 1 ) = SQRT( ( MLEWD( 2 , 1) - ( XNREP*( MLEWD( 3 , 1)**2) ) ) / ( XNREP - 1 ) ) 
MLESP( 4 , I ) = SQRT( (MLESP( 2 , I ) -(XNREP*(MLESP( 3 , I )**2) ) ) /(XNREP- 1) ) 
REGESTC 4 , 1 )=SQRT( ( REGEST( 2 , 1 ) - ( XNREP*( REGEST( 3 , 1 )**2 ) ) ) / ( XNREP- 1 ) ) 
WRESTK 4 , 1 )=SQRT( (WRESTK 2 , 1 ) -(XNREP-'-( WRESTK 3 , 1 )**2) ) )/(XNREP- 1) ) 
WREST2C 4 , 1 )=SQRT( ( WREST2 ( 2 , 1 ) - ( XNREP*( WREST2( 3 , 1 )**2 ) ) ) / ( XNREP- 1) ) 
V/REST3( 4 , 1 )=SQRT( ( WREST3( 2 , 1 ) - ( XNREP*( WREST3( 3 , 1 )**2 ) ) ) / ( XNREP - 1 ) ) 
600 CONTINUE 



C PRINT THE FINAL OUTPUT TABLE TO A FILE 



601 WRITE(30,3000) 

3000 FORMATC ' O' ,T47,' DISCRETE RELIABILITY GROWTH SIMULATION') 
WRITE(30,3010) 

3010 FORMAT( '-' ,T54, 'MODEL PARAMETER SUMMARY') 

WRITE( 30,3020) NCAUSE 

3020 FORMATC 'o' ,T4 7, 'NUMBER OF POSSIBLE FAILURE CAUSES ',14) 

IF (FRELOP .EQ. 1) GOTO 4000 
WRITE(30,3030) 

3030 FORMATC 'O' ,T38, 'CAUSE NUMBER' ,T64 ,' SINGLE TRIAL PR( SUCCESS ) FOR 
CPHASE 1') 

DO 3050 M=l, NCAUSE 
WRITE(30,3040) M,A(M,1) 

3040 FORMATC' ' ,T43 , 12 ,T79 ,F8. 6) 

3050 CONTINUE 

WRITE(30,3060) FRIMP 

3060 FORMATC 'o' ,T37, 'FRACTION CAUSE RELIABILITY IMPROVES AFTER FAILURE 
C ,F8. 6) 

5000 WRITEC 30,3080) NPHASE 
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3080 FORMATC ,T48,' NUMBER OF PHASES IN THE SIMULATION ’,12) 
WRITE(30,3090) 

3090 FORMATC ' O' ,T42, 'PHASE NUMBER' ,T59 ,' NUMBER OF FAILURES IN THE FIRST 
C PHASE') 

DO 3110 M=1,NPHASE 
WRITE(30,3100) M,NFAPH(M) 

3100 FORMATC' ' ,T43 , I2,T73 , 12) 

3110 CONTINUE 

WRITE(30,3120) NFAIL 

3120 FORMATC 'O' ,T51, 'TOTAL NUMBER OF FAILURES ',14) 

IFCDISOPT. EQ. 2) GO TO 3160 
WRITE(30,3130) 

3130 FORMATC '-' ,T38, 'DISCOUNTING PERFORMED USING THE CONSTANT FRACTION 
CMETHOD' ) 

WRITEC30,3140) R 

3140 FORMATC 'O' ,T44, 'FRACTION EACH FAILURE IS DISCOUNTED ' ,F8. 6) 
WRITEC30,3150) N 

3150 FORMATC' ' ,T33 ,' NUMBER OF TRIALS AFTER A FAILURE BEFORE A DISCOUNT 
C IS APPLIED ' ,14) 

GO TO 3190 
3160 WRITEC30,3170) 

3170 FORMATC '-' ,T44, 'DISCOUNTING PERFORMED USING THE LLOYD METHOD') 
WRITEC30,3180) GAMA 

3180 FORMATC 'O' ,T39, 'PERCENT C. I. C USED AS DISCOUNT FRACTION ) ' ,F8. 6 
C) 

WRITEC30,3185) LDI 

3185 FORMATC 'o' ,T50, 'LLOYD DISCOUNT INTERVAL: ',13) 

3190 WRITEC30,3200) DSEEDl 

3200 FORMATC '-' ,T46, 'RANDOM NUMBER SEED USED ' ,F15. 2) 

WRITEC30,3210) NREP 

3210 FORMATC 'O' ,T37, 'NUMBER OF REPETITIONS OF THE SIMULATION PERFORMED 
C ,17) 

WRITEC30,3220) 

3220 FORMATC ' l' ,T61, 'ESTIMATOR: ' ) 

WRITEC30,3230) 

3230 FORMATC 'O' ,T48, 'SINGLE PHASE MLE WITHOUT DISCOUNTING') 
WRITEC30,3240) 

3240 FORMATC '-' ,T60, 'MEAN' ,T83, 'ESTIMATE' ,T109, '95 %' ) 

WRITEC30,3250) 

3250 FORMATC' ',T12, 'PHASE NUMBER' ,T29 ,' ACTUAL RELIABILITY' ,T52 ,' PREDIC 
CTED RELIABILITY' ,T78,' STANDARD DEVIATION' ,T101 ,' CONFIDENCE INTERVA 
CL' ) 

C COMPUTE C. I. FOR SINGLE PHASE MLE 

DO 3270 M=1,NPHASE 

Cl = Cl. 96*MLESPC4,M))/SQRTCXNREP) 

CIU = MLESPC3,M) + Cl 
CIL = MLESPC3,M) - Cl 

WRITEC30,3260) M, ARELCM) ,MLESPC 3 ,M) ,MLESPC4 ,M) ,CIL,CIU 
3260 FORMATC 'o' ,T17,I2,T34,F8. 6,T58,F8. 6,T82,F9. 6,T99, 'C ',F8.6,' , ' ,F 
C8.6,' )') 

3270 CONTINUE 

WRITEC 30,3220) 

WRITEC30,3280) 

3280 FORMATC ' o' ,T42, 'MAX LIKELIHOOD ESTIMATE USING DISCOUNTED FAILURES' 
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C) 

WRITE(30,3240) 

WRITE(30,3250) 

C COMPUTE C. I. FOR MLE WITH DISCOUNTING 



DO 3290 M=1,NPHASE 

Cl = (1. 96*MLEWD(4,M))/SQRT(XNREP) 

CIU = MLEWD(3,M) + Cl 
CIL = MLEWD(3,M) - Cl 

WT?ITE(30,3260) M,AREL(M) ,MLEWD( 3 ,M) ,MLEV,T)( 4 ,M) , CIL, CIU 
3290 CONTINUE 

WRITE(30,3220) 

WRITE(30,3300) 

3300 FORMATC’O’ ,T3 8, 'WEIGHTED AVERAGE ESTIMATE USING FAILURE DISCOUNTIN 
CG’ ) 

WRITE(30,3240) 

WRITE(30,3250) 

C COMPUTE C. I. FOR WOODS WEIGHTED AVERAGE ESTIMATES 



DO 3310 M=1,NPHASE 

Cl = (1. 96*EST(4,M))/SQRT(XNT1EP) 

CIU = EST(3,M) + Cl 
CIL = EST(3,M) - Cl 

WRITE(30,3260) M, AREL(M) ,EST( 3 ,M) ,EST(4,M) , CIL, CIU 
3310 CONTINUE 



WRITE(30,3220) 

WRITE(30,3320) 

3320 FORMATC 'O' ,T43,’ REGRESSION ESTIMATE USING DISCOUNTED FAILURES’) 
WRITE(30,3240) 

WRITE(30,3250) 

C COMPUTE C.I. FOR EXPONENTIAL REGRESSION ESTIMATES 



DO 3330 M=1,NPHASE 

Cl = (1. 96*REGEST(4,M))/SQRT(XNREP) 

CIU = REGEST(3,M) + Cl 
CIL = REGEST(3,M) - Cl 

WRITE(30,3260) M,AREL(M) ,REGEST(3 ,M) ,REGEST(4,M) , CIL, CIU 
3330 CONTINUE 



C WEIGHTED REGRESSION (METHOD 1) 

WRITE(30,3220) 

WRITE(30,3321) 

3321 FORMATC 'O' ,T4 3, 'WEIGHTED REGRESSSION ESTIMATE (METHOD 1 )') 

WRITE(30,3240) 

WRITE(30,3250) 

C COMPUTE C. I. FOR WEIGHTED EXPONENTIAL REGRESSION ESTIMATES 



DO 3331 M=1,NPHASE 

Cl = (1. 96*WREST1(4,M))/SQRT(XNREP) 
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CIU = WREST1(3,M) + Cl 
CIL = WREST1(3,M) - Cl 

WRITE(30,3260) M,AREL(M) ,WREST1(3,M) ,WREST1(4,M) , CIL, CIU 

3331 CONTINUE 

C WEIGHTED REGRESSION (METHOD 2) 

WRITE(30,3220) 

WRITE(30,3322) 

3322 FORMATC 'O' ,T4 3, 'WEIGHTED REGRESSION ESTIMATE (METHOD 2) ') 

WRITE(30,3240) 

WRITE(30,3250) 

C COMPUTE C. I. FOR WEIGHTED EXPONENTIAL REGRESSION ESTIMATES 

DO 3332 M=1,NPHASE 

Cl = (1. 96*WREST2(4,M))/SQRT(XNREP) 

CIU = WREST2(3,M) + Cl 
CIL = WREST2(3,M) - Cl 

WRITE(30,3260) M,AREL(M) ,WREST2(3,M) ,WREST2(4,M) ,CIL,CIU 

3332 CONTINUE 

C WEIGHTED REGRESSION (METHOD 3) 

WRITE(30,3220) 

WRITE(30,3323) 

3323 FORMAT( 'O' ,T43, 'WEIGHTED REGRESSION ESTIMATE (METHOD 3) ') 

WRITE(30,3240) 

WRITE(30,3250) 

C COMPUTE C. I. FOR WEIGHTED EXPONENTIAL REGRESSION ESTIMATES 



DO 3333 M=1,NPHASE 

Cl = (1. 96*WREST3(4,M))/SQRT(XNREP) 

CIU = WREST3(3,M) + Cl 
CIL = WREST3(3,M) - Cl 

WRITE(30,3260) M,AREL(M) ,WREST3( 3 ,M) ,WREST3(4,M) , CIL, CIU 
3333 CONTINUE 

WRITE(30,3340) 

3340 FORMATC 1' ,T59 ,' RECAPITULATION' // ) 

WRITE(30,3350) 

3350 FORMAT( ' -' ,T3, 'PHASE' , Til, 'ACTUAL' ,T28, 'MEAN' ,T38, 'EST' ,T50, 'MEAN' 
C,T62,'EST' ,T72,'MEAN' ,T82,'EST' ,T92, 'MEAN' ,T102,'EST' ,T112 
C,'MEAN' ,T122,'EST') 

WRITE(30,3360) 

3360 FORMATC ' ,T11 ,'RELIAB' ,T28,'WGT' ,T38,'STD' ,T50,'MLE' ,T62,'STD' ,T7 
C2, 'PHASE' ,T82, ^STD' ,T92, 'REG' ,T102, 'STD' ,T112, 'WT' ,T122,'STD' ) 
WRITE(30,3370) 

3370 FORMATC ' ,T28 , ' AVG' ,T38 , ' DEV' ,T50 , ' W/D ' ,T62 , ' DEV' ,T72 

C, 'MLE' ,T82, 'DEV' ,T92, 'EST' ,T102, 'DEV' ,T112, 'REG' ,T122, 'DEV' ) 
WRITE(30,3375) 

3375 FORMATC ' ,T28 , ' EST' ,T112 , ' ( 1) ' / ) 

DO 650 I=1,NPHASE 

WRITE(30,3380)I,AREL(I) ,EST(3,I) ,EST(4,I) ,MLEWD(3,I) ,MLEWD(4,I) , 
CMLESP( 3,1), MLESP(4 , I ) , REGEST( 3,1), REGEST( 4,1) , WREST1( 3 , I ) , 
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CWREST1(4,I) 

3380 FORMATC 'O' ,T4 , 13 ,T11 ,F7. 6 ,T26 ,F7. 6 ,T36 ,F7. 6 ,T50,F7. 6,T62,F7. 6,T72, 
CF7. 6,T82,F7. 6,T92,F7. 6,T102,F7. 6,T112,F7. 6,T122,F7. 6) 

650 CONTINUE 

C PAGE 2 OF RECAPITULATION: 

WRITE(30,3381) 

3381 FORMATC ' 1' ,T59, 'RECAPITULATION CONT. '//) 

WRITE(30,3382) 

3382 FORMATC ' ,T3 ,' PHASE ' .Til ,' ACTUAL' ,T28,'MEAN' ,T38, 'EST' T53, 
C'MEAN' ,T63,'EST' ,T78/MEAN' ,T88,'EST' ,T103,'MEAN' ,T113/EST') 

WRITEC30,3383) 

3383 FORMATC’ ' ,Tli , ' RELIAB' .T28.'REG' ,T38,’STD' , T5 3 ,' METHOD ' ,T63, 
C'STD' ,T77/METHOD' ,T88, 'STD^ ,T103, 'METHOD' ,T113, 'STD' ) 

WRITEC30,3384) 

3384 FORMATC ' ' ,T28, 'EST' ,T35 ,' DEVIATION' ,T53, 'ONE' ,T60 ,' DEVIATION' 
C,T78, 'TWO' ,T85, 'DEVIATION' ,T103, 'THREE' ,T110, 'DEVIATION' ) 

i)0 651 I=1,NPHASE 

WRITEC30.3385) I , ARELC I ) ,REGESTC 3 , I) ,REGESTC4 , I) ,WREST1C 3 , I) 

C , WRESTIC 4,1) , WREST2C 3,1) , WREST2C 4 , 1 ) ,WREST3C 3,1) ,WREST3C 4 , 1 ) 

3385 FORMATC 'O' ,T4 , 13 ,T11 ,F7. 6 ,T26 ,F7. 6 ,T36 ,F7. 6 ,T51 ,F7. 6 ,T61 ,F7. 6 
C,T76,F7. 6,T86,F7. 6,T101,F7. 6,T111,F7. 6) 

WRITE C 5 3, 3400) MLEWDC4,!) 

WRITE C 54, 3400) MLEWDC3,I) 

WRITE C 63, 3400) MLESPC4,!) 

WRITE C 64, 3400) MLESPC3,!) 

WRITE C 73, 3400) REGESTC4,!) 

WRITE C 74, 3400) REGESTC3,I) 

WRITEC 17,3400) WRESTIC 4, I) 

WRITEC 18,3400) WRESTIC 3, I) 

WRITEC 9 3, 3400) WREST2C4,!) 

WRITEC 94, 3400) WREST2C3,!) 

WRITEC 95, 3400) WREST3C4,!) 

WRITEC 96, 3400) WREST3C3,!) 

WRITE C 84, 3400) ARELC I) 

3400 FORMATC ' ’ ,F7. 6) 

651 CONTINUE 



GO TO 6000 
4000 WRITEC 30,4010) 

4010 FORMATC IX, //,T50, 'FIXED PHASE RELIABILITY OPTION') 
WRITEC 30,4020) 

4020 FORMATC ' -' ,T38, 'PHASE NUMBER' ,T78, 'ACTUAL RELIABILITY' ) 
DO 4030 M=1,NPHASE 
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WRITE(30,4040) M,AREL(M) 
4040 FORMATC 'O' ,T41 , 12 ,T83 ,F8. 6) 
4030 CONTINUE 
GO TO 5000 
6000 CONTINUE 

STOP 

END 
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5. SAMPLE OUTPUT FROM DRG FORTRAN A1 

DISCRETE RELIABILITY GROWTH SIMULATION 

MODEL PARAMETER SUMMARY 
NUMBER OF POSSIBLE FAILURE CAUSES 5 

FIXED PHASE RELIABILITY OPTION 



PHASE NUMBER 
1 
2 

3 

4 

5 

6 

7 

8 
9 

10 



ACTUAL RELIABILITY 
0.398418 
0.428109 
0.480793 
0.539243 
0.609949 
0.702683 
0.798124 
0.899963 
0.950990 
0.990040 



PHASE 

1 

2 

3 

4 

5 

6 

7 

8 
9 

10 



NUMBER OF PHASES IN THE SIMULATION 10 

NUMBER NUMBER OF FAILURES IN THE FIRST PHASE 

1 

1 

1 

1 

1 

1 

1 

1 

1 

1 

TOTAL NUMBER OF FAILURES 10 



DISCOUNTING PERFORMED USING THE CONSTANT FRACTION METHOD 

FRACTION EACH FAILURE IS DISCOUNTED 0.000000 
NUMBER OF TRIALS AFTER A FAILURE BEFORE A DISCOUNT IS APPLIED 1 



RANDOM NUMBER SEED USED 624712.00 

NUMBER OF REPETITIONS OF THE SIMULATION PERFORMED 500 
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ESTIMATOR t 

MAX LIKELIHOOD ESTIMATE USING DISCOUNTED FAILURES 
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ESTIMATOR! 

REGRESSION ESTIMATE USING DISCOUNTED FAILURES 
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